我要投稿
投稿指南
RSS订阅
1.1.3. 定义一个投影系统
一个投影坐标系统（如UTM等）基于一个地理坐标系统定义以便于在线性单位和角度位置之间转换。接下来的代码定义了一个在WGS84基准面下的地理坐标系统为基础的UTM17带投影坐标系统。
1 OGRSpatialReference oSRS;
2 oSRS.SetProjCS( "UTM 17 (WGS84) in northern hemisphere." );
3 oSRS.SetWellKnownGeogCS( "WGS84" );
4 oSRS.SetUTM( 17, TRUE );
调用SetProjCS()设置一个用户定义名字的投影坐标系统并确定系统被投影过。SetWellKnownGeogCS()分配一个地理坐标系统，SetUTM()设置投影转换参数细节。上面的步骤直到这时才能创建一个合法的定义。但是这个对象需要时将会重新自动编排内在表达，以保持合法性。
请尤其注意定义对象时的步骤顺序！！
上面的定义会定义一个像下面一样的WKT。注意UTM117被扩展成横轴莫卡托定义UTM带的参数。
1 PROJCS["UTM 17 (WGS84) in northern hemisphere.",
2 GEOGCS["WGS 84",
3 DATUM["WGS_1984",
4 SPHEROID["WGS 84",6378137,298.257223563,
5 AUTHORITY["EPSG",7030]],
6 TOWGS84[0,0,0,0,0,0,0],
7 AUTHORITY["EPSG",6326]],
8 PRIMEM["Greenwich",0,AUTHORITY["EPSG",8901]],
9 UNIT["DMSH",0.0174532925199433,AUTHORITY["EPSG",9108]],
10 AXIS["Lat",NORTH],
11 AXIS["Long",EAST],
12 AUTHORITY["EPSG",4326]],
13 PROJECTION["Transverse_Mercator"],
14 PARAMETER["latitude_of_origin",0],
15 PARAMETER["central_meridian",-81],
16 PARAMETER["scale_factor",0.9996],
17 PARAMETER["false_easting",500000],
18 PARAMETER["false_northing",0]]
不同的投影会有不同的函数来设定参数(SetTM() (Transverse Mercator), SetLCC() (Lambert Conformal Conic)和 SetMercator())
1.1.4. 查询坐标系统
一旦OGRSpatialReference建立起来了，就可以查询多种信息。用OGRSpatialReference::IsProjected() 和 OGRSpatialReference::IsGeographic()可以看他是否是投影坐标还是地理坐标。OGRSpatialReference::GetSemiMajor(), OGRSpatialReference::GetSemiMinor() 和OGRSpatialReference::GetInvFlattening()可以获取椭球体信息。OGRSpatialReference::GetAttrValue() 可以用来获取PROJCS, GEOGCS, DATUM, SPHEROID, 和PROJECTION 的名字表达字符串。OGRSpatialReference::GetProjParm()可以用来获取投影参数。OGRSpatialReference::GetLinearUnits() 可以用来获取线性单位，并转换成米。
下面的代码示范了GetAttrValue()获取投影的用法，以及GetProjParm()GetAttrValue() 的用法。已经定义的投影参数(如SRS_PP_CENTRAL_MERIDIAN)应该在GetProjParm()获取投影参数的时候使用。ogrspatialreference.cpp中设置不同投影的方法的代码可以用来查找哪个投影用哪个参数。
1 const char *pszProjection = poSRS->GetAttrValue("PROJECTION");
2 if( pszProjection == NULL )
3 {
4 if( poSRS->IsGeographic() )
5 sprintf( szProj4+strlen(szProj4), "+proj=longlat " );
6 else
7 sprintf( szProj4+strlen(szProj4), "unknown " );
8 }
9 else if( EQUAL(pszProjection,SRS_PT_CYLINDRICAL_EQUAL_AREA) )
10 {
11 sprintf( szProj4+strlen(szProj4),
12 "+proj=cea +lon_0=%.9f +lat_ts=%.9f +x_0=%.3f +y_0=%.3f ",
13 poSRS->GetProjParm(SRS_PP_CENTRAL_MERIDIAN,0.0),
14 poSRS->GetProjParm(SRS_PP_STANDARD_PARALLEL_1,0.0),
15 poSRS->GetProjParm(SRS_PP_FALSE_EASTING,0.0),
16 poSRS->GetProjParm(SRS_PP_FALSE_NORTHING,0.0) );
17 }
18 ...
1.1.5. 坐标转换
OGRCoordinateTransformation 类用来在不同坐标系统间转换坐标。新的转换对象由OGRCoordinateTransformation()创建。创建后OGRCoordinateTransformation::Transform() 方法就可以用来转换不同坐标系的点了。
1 OGRSpatialReference oSourceSRS, oTargetSRS;
2 OGRCoordinateTransformation *poCT;
3 double x, y;
4 oSourceSRS.importFromEPSG( atoi(papszArgv[i+1]) );
5 oTargetSRS.importFromEPSG( atoi(papszArgv[i+2]) );
6 poCT = OGRCreateCoordinateTransformation( &oSourceSRS,
7 &oTargetSRS );
8 x = atof( papszArgv[i+3] );
9 y = atof( papszArgv[i+4] );
10 if( poCT == NULL || !poCT->Transform( 1, &x, &y ) )
11 printf( "Transformation failed.\n" );
12 else
13 printf( "(%f,%f) -> (%f,%f)\n",
14 atof( papszArgv[i+3] ),
15 atof( papszArgv[i+4] ),
16 x, y );
有两个原因导致在转换时会出错。第一，OGRCreateCoordinateTransformation() 可能失败。一般是因为程序认为两个指定坐标系统不能建立转换关系，这可能是由于使用的投影Proj内部暂时不支持，这样两个不一样的基准面就没有已知关系。或者一个投影定义得不适当。如果OGRCreateCoordinateTransformation() 失败了将返回空。
OGRCoordinateTransformation::Transform() 方法本身也可能失败。这可能是因为一个有上面问题积累后形成的问题。或者是一个“因为一个或者多个点因数值上未定义而省略操作”起因的后果。Transform()函数如果成功，返回TRUE，如果某个点转换失败，剩下的点就会以一个错误的状态保持。
坐标转换服务其实可以处理3D点的。并且服务会根据椭球体和基准面上的高程差异调节高程。将来，也有可能会有在不同矢量基准面上的点间转换的应用。如果没有Z值，则假设点都在大地水准面上。
下面的例子显示了如何方便地创建一个地理坐标系统到一个基于同一个地理坐标系统的投影坐标系统。并在两个坐标系统之间转换。
1 OGRSpatialReference oUTM, *poLatLong;
2 OGRCoordinateTransformation *poTransform;
3 oUTM.SetProjCS("UTM 17 / WGS84");
4 oUTM.SetWellKnownGeogCS( "WGS84" );
5 oUTM.SetUTM( 17 );
6 poLatLong = oUTM.CloneGeogCS();
7 poTransform = OGRCreateCoordinateTransformation( &oUTM, poLatLong );
8 if( poTransform == NULL )
9 {
10 ...
11 }
12 ...
13 if( !poTransform->Transform( nPoints, x, y, z ) )
14 ...
1.1.6. API绑定
C的接口在ogr_srs_api.h中定义，python绑定在osr.py模块中。一些C++的方法在C和python绑定中没有定义。
Python Bindings(C绑定略过)
1 class osr.SpatialReference
2 def __init__(self,obj=None):
3 def ImportFromWkt( self, wkt ):
4 def ExportToWkt(self):
5 def ImportFromEPSG(self,code):
6 def IsGeographic(self):
7 def IsProjected(self):
8 def GetAttrValue(self, name, child = 0):
9 def SetAttrValue(self, name, value):
10 def SetWellKnownGeogCS(self, name):
11 def SetProjCS(self, name = "unnamed" ):
12 def IsSameGeogCS(self, other):
13 def IsSame(self, other):
14 def SetLinearUnits(self, units_name, to_meters ):
15 def SetUTM(self, zone, is_north = 1):
16 class CoordinateTransformation:
17 def __init__(self,source,target):
18 def TransformPoint(self, x, y, z = 0):
19 def TransformPoints(self, points):
1.1.7. 接口实现
OGRCoordinateTransformation服务的实现是建立在PROJ4库之上，该库是最先由USGS的Gerald Evenden 编写的。
