GDAL线面互转换(2)
在上一个文章中介绍了线转化为面和面转化为线,其主要的实现思路就是把面中的点取出来构成线,把线中的点取出来构成面,实际上就是一个硬拷贝,无奈客户的实际需求并非如此,客户想要线转面的时候几条相交线构成面,面转线的时候相同的线去除,所以又重新对功能进行了调整。
关于线构面,这个过程中也是挺曲折的,在GDAL的群里问了好久,大家给的答案一致是需要自己写,线构面的算法需要使用左转或者右转算法,并且在网上查了相关的资料,发现蒋波涛编著的《插件式GIS应用框架的设计与实现:基于C#和AE 9.2》一书中有ArcGIS实现GDAL的完整代码,没办法,埋头去干吧,先做了一部分前期的验证,验证的时候,公司一大牛问我OGRGeometry中的Polygonize ()方法是干啥的?仔细查阅相关资料,就是要做线转面的啊!赶紧验证,发现失败,后来又仔细分析,似乎该方法需要OGRGeometry的空间拓扑关系正确,如何保证这点儿呢?考虑用Union来把一根根线合并了拓扑关系应该没有问题吧,实际测试了一下,的确正确,心中窃喜。后来拿实际的省界线来合成省面,发现效率超级慢,跟踪代码发现,OGRGeometry合并的时候随着合并线段的增加时间也在快速的增加……问题找到,上网查找解决方案,看到两篇文章(http://www.docin.com/p-1155026547.html)(http://www.docin.com/p-1238395230.html)介绍说可以考虑使用UnionCascaded(级联求并)可以大大的加快效率,实际验证确不知道怎么使用,考虑到项目的进度要求,这个问题暂且搁置了,后续再提高吧。
对于面构线,就相对简单了不少,先把一个个的面转化成线(硬拷贝,参考上一博文中的方法),再把转化成的线Union了,生成一个拓扑关系正确的OGRMultiPolygon,再调用Simplify方法,得到的结果即为想要的线,不过当线特别多还很复杂的时候效率也高,问题和线构面的时候一样。
在此补充上源代码,在此留个记录,后续考虑怎么解决这个问题吧……
![](https://images.cnblogs.com/OutliningIndicators/ContractedBlock.gif)
1 /* 2 * @brief ConvertPolygonToPolyline 面图层转换为线图层 3 * @param[in] tString polylinePath 转换后线图层文件路径 4 * @param[in] OGRLayer* pLayer 要转换的面图层文件 5 * @param[in] Envelope envelope 要转换数据的范围 6 * @param[in] vector<long> vecFIDs 选中的要素ID列表 7 * @param[in] pOGRSpatialReference 要转换的数据的空间参考(如果为空表示坐标系信息不变) 8 * @return bool 是否成功 9 * @author 10 * @date 11 * @note 2015年11月04日 小八创建; 12 */ 13 bool FeatureLayerOperator::ConvertPolygonToPolylineEx(tString polylinePath,OGRLayer* pLayer,Envelope envelope,vector<long> vecFIDs,OGRSpatialReference* pOGRSpatialReference) 14 { 15 // 判断 16 if(pLayer==NULL) return false; 17 if(true==polylinePath.empty()) return false; 18 19 // 坐标系读取 20 OGRSpatialReference* pOGRSpatialReference_Source=pLayer->GetSpatialRef(); 21 bool isSameCoordSystem = false; 22 if(pOGRSpatialReference == NULL) 23 { 24 pOGRSpatialReference=pOGRSpatialReference_Source; 25 isSameCoordSystem=true; 26 } 27 else if(pOGRSpatialReference!=NULL && pOGRSpatialReference_Source!=NULL) 28 { 29 isSameCoordSystem=pOGRSpatialReference_Source->IsSame(pOGRSpatialReference); 30 } 31 32 // 创建Shape文件 33 OGRDataSource* pOGRDataSource=CreateShapeFile(polylinePath,pOGRSpatialReference,wkbLineString); 34 if(pOGRDataSource==NULL) return false; 35 36 OGRLayer* pOGRLayer=pOGRDataSource->GetLayer(0); 37 if(pOGRLayer==NULL) return false; 38 39 // 面转线再合并 40 OGRFeature* pOGRFeature_Old; 41 OGRGeometry* pTempGeometry=NULL; 42 OGRGeometry* pTempGeometryUnion=NULL; 43 44 // 当前选择导出 45 if(false==vecFIDs.empty()&&vecFIDs.size()>0) 46 { 47 for(int i=0;i<vecFIDs.size();i++) 48 { 49 pOGRFeature_Old=pLayer->GetFeature(vecFIDs[i]); 50 pTempGeometry=ConvertPolygonToPolylineGeo(pOGRFeature_Old->GetGeometryRef()); 51 pTempGeometry->assignSpatialReference(pOGRSpatialReference_Source); 52 if(false == isSameCoordSystem)pTempGeometry->transformTo(pOGRSpatialReference); 53 54 if(pTempGeometryUnion==NULL) pTempGeometryUnion=pTempGeometry; 55 else pTempGeometryUnion=pTempGeometryUnion->Union(pTempGeometry); 56 } 57 } 58 else 59 { 60 if(false!=envelope.isNull()&&envelope.getMaxX()!=envelope.getMinX()) 61 { 62 pLayer->SetSpatialFilterRect(envelope.getMinX(),envelope.getMinY(),envelope.getMaxX(),envelope.getMaxY()); 63 } 64 65 pOGRFeature_Old=pLayer->GetNextFeature(); 66 while(NULL!= pOGRFeature_Old) 67 { 68 pTempGeometry=ConvertPolygonToPolylineGeo(pOGRFeature_Old->GetGeometryRef()); 69 pTempGeometry->assignSpatialReference(pOGRSpatialReference_Source); 70 if(false == isSameCoordSystem)pTempGeometry->transformTo(pOGRSpatialReference); 71 72 if(pTempGeometryUnion==NULL) pTempGeometryUnion=pTempGeometry; 73 else pTempGeometryUnion=pTempGeometryUnion->Union(pTempGeometry); 74 75 pOGRFeature_Old=pLayer->GetNextFeature(); 76 } 77 } 78 79 // 获得Simply的距离 80 double distanceValue=0.0; 81 if(NULL==pOGRSpatialReference) // 如果为空 82 { 83 OGREnvelope pTempOGREnvelope ; 84 pTempGeometryUnion->getEnvelope(&pTempOGREnvelope); 85 if(pTempOGREnvelope.MaxX<180) distanceValue=0.00000001; 86 else distanceValue=0.01; 87 } 88 else if(true==pOGRSpatialReference->IsProjected()) // 如果是Project的 89 { 90 distanceValue=0.01; 91 } 92 else // 如果是Geo的 93 { 94 distanceValue=0.00000001; 95 } 96 pTempGeometryUnion=pTempGeometryUnion->Simplify(distanceValue); 97 98 // 在ShapeFile文件中添加数据行 99 OGRFeature* pOGRFeature_New; 100 OGRGeometry* pOGRGeometry; 101 OGRFeatureDefn* pOGRFeatureDefn=NULL; 102 pOGRFeatureDefn=pOGRLayer->GetLayerDefn(); 103 104 OGRwkbGeometryType ogrGeometryType=pTempGeometryUnion->getGeometryType(); 105 ogrGeometryType=wkbFlatten(ogrGeometryType); 106 107 if(ogrGeometryType==OGRwkbGeometryType::wkbMultiLineString) 108 { 109 OGRGeometryCollection* pOGRGeometryCollectionTarget=(OGRGeometryCollection*) pTempGeometryUnion; 110 int geometryCount=pOGRGeometryCollectionTarget->getNumGeometries(); 111 for(int i=0;i<geometryCount;i++) 112 { 113 pOGRFeature_New=OGRFeature::CreateFeature(pOGRFeatureDefn); 114 pOGRGeometry=pOGRGeometryCollectionTarget->getGeometryRef(i); 115 pOGRFeature_New->SetGeometry(pOGRGeometry); 116 pOGRLayer->CreateFeature(pOGRFeature_New); 117 118 OGRFeature::DestroyFeature(pOGRFeature_New); 119 pOGRFeature_New=NULL; 120 } 121 } 122 else if(ogrGeometryType==OGRwkbGeometryType::wkbLineString) 123 { 124 pOGRFeature_New=OGRFeature::CreateFeature(pOGRFeatureDefn); 125 pOGRGeometry=pTempGeometryUnion; 126 pOGRFeature_New->SetGeometry(pOGRGeometry); 127 pOGRLayer->CreateFeature(pOGRFeature_New); 128 129 OGRFeature::DestroyFeature(pOGRFeature_New); 130 pOGRFeature_New=NULL; 131 132 } 133 OGRDataSource::DestroyDataSource(pOGRDataSource); 134 135 // 销毁pTargetGeometrys 136 OGRGeometryFactory::destroyGeometry(pTempGeometryUnion); 137 pTempGeometryUnion=NULL; 138 139 return true; 140 } 141 142 /* 143 * @brief ConvertPolygonToPolyline 线图层转换为面图层 144 * @param[in] tString polylinePath 转换后面图层文件路径 145 * @param[in] OGRLayer* pLayer 要转换的线图层文件 146 * @param[in] Envelope envelope 要转换数据的范围 147 * @param[in] vector<long> vecFIDs 选中的要素ID列表 148 * @param[in] pOGRSpatialReference 要转换的数据的空间参考(如果为空表示坐标系信息不变) 149 * @return bool 是否成功 150 * @author 151 * @date 152 * @note 2015年11月04日 小八创建; 153 */ 154 bool FeatureLayerOperator::ConvertPolylineToPolygonEx(tString polylinePath,OGRLayer* pLayer,Envelope envelope,vector<long> vecFIDs,OGRSpatialReference* pOGRSpatialReference) 155 { 156 // 判断 157 if(pLayer==NULL) return false; 158 if(true==polylinePath.empty()) return false; 159 160 // 坐标系读取 161 OGRSpatialReference* pOGRSpatialReference_Source=pLayer->GetSpatialRef(); 162 bool isSameCoordSystem = false; 163 if(pOGRSpatialReference == NULL) 164 { 165 pOGRSpatialReference=pOGRSpatialReference_Source; 166 isSameCoordSystem=true; 167 } 168 else if(pOGRSpatialReference!=NULL && pOGRSpatialReference_Source!=NULL) 169 { 170 isSameCoordSystem=pOGRSpatialReference_Source->IsSame(pOGRSpatialReference); 171 } 172 173 // 创建Shape文件 174 OGRDataSource* pOGRDataSource=CreateShapeFile(polylinePath,pOGRSpatialReference,wkbPolygon); 175 if(pOGRDataSource==NULL) return false; 176 177 OGRLayer* pOGRLayer=pOGRDataSource->GetLayer(0); 178 if(pOGRLayer==NULL) return false; 179 180 // 面合并 181 OGRFeature* pOGRFeature_Old; 182 OGRGeometry* pTempGeometry=NULL; 183 OGRGeometry* pTempGeometryUnion=NULL; 184 185 // 当前选择导出 186 if(false==vecFIDs.empty()&&vecFIDs.size()>0) 187 { 188 for(int i=0;i<vecFIDs.size();i++) 189 { 190 pOGRFeature_Old=pLayer->GetFeature(vecFIDs[i]); 191 pTempGeometry=pOGRFeature_Old->GetGeometryRef(); 192 if(false == isSameCoordSystem)pTempGeometry->transformTo(pOGRSpatialReference); 193 194 if(pTempGeometryUnion==NULL) pTempGeometryUnion=pTempGeometry; 195 else pTempGeometryUnion=pTempGeometryUnion->Union(pTempGeometry); 196 } 197 } 198 else 199 { 200 if(false!=envelope.isNull()&&envelope.getMaxX()!=envelope.getMinX()) 201 { 202 pLayer->SetSpatialFilterRect(envelope.getMinX(),envelope.getMinY(),envelope.getMaxX(),envelope.getMaxY()); 203 } 204 205 pOGRFeature_Old=pLayer->GetNextFeature(); 206 while(NULL!= pOGRFeature_Old) 207 { 208 pTempGeometry=pOGRFeature_Old->GetGeometryRef(); 209 if(false == isSameCoordSystem)pTempGeometry->transformTo(pOGRSpatialReference); 210 211 if(pTempGeometryUnion==NULL) pTempGeometryUnion=pTempGeometry; 212 else pTempGeometryUnion=pTempGeometryUnion->Union(pTempGeometry); 213 214 pOGRFeature_Old=pLayer->GetNextFeature(); 215 } 216 } 217 218 OGRGeometry* pOGRGeometryUnion = pTempGeometryUnion->Polygonize(); 219 220 // 在ShapeFile文件中添加数据行 221 OGRFeature* pOGRFeature_New; 222 OGRGeometry* pOGRGeometry; 223 OGRFeatureDefn* pOGRFeatureDefn=NULL; 224 pOGRFeatureDefn=pOGRLayer->GetLayerDefn(); 225 226 OGRwkbGeometryType ogrGeometryType=pOGRGeometryUnion->getGeometryType(); 227 ogrGeometryType=wkbFlatten(ogrGeometryType); 228 229 if(ogrGeometryType==OGRwkbGeometryType::wkbGeometryCollection||ogrGeometryType==OGRwkbGeometryType::wkbMultiPolygon) 230 { 231 OGRGeometryCollection* pOGRGeometryCollectionTarget=(OGRGeometryCollection*) pOGRGeometryUnion; 232 int geometryCount=pOGRGeometryCollectionTarget->getNumGeometries(); 233 for(int i=0;i<geometryCount;i++) 234 { 235 pOGRFeature_New=OGRFeature::CreateFeature(pOGRFeatureDefn); 236 pOGRGeometry=pOGRGeometryCollectionTarget->getGeometryRef(i); 237 pOGRFeature_New->SetGeometry(pOGRGeometry); 238 pOGRLayer->CreateFeature(pOGRFeature_New); 239 240 OGRFeature::DestroyFeature(pOGRFeature_New); 241 pOGRFeature_New=NULL; 242 } 243 } 244 else if(ogrGeometryType==OGRwkbGeometryType::wkbPolygon) 245 { 246 pOGRFeature_New=OGRFeature::CreateFeature(pOGRFeatureDefn); 247 pOGRGeometry=pOGRGeometryUnion; 248 pOGRFeature_New->SetGeometry(pOGRGeometry); 249 pOGRLayer->CreateFeature(pOGRFeature_New); 250 251 OGRFeature::DestroyFeature(pOGRFeature_New); 252 pOGRFeature_New=NULL; 253 254 } 255 OGRDataSource::DestroyDataSource(pOGRDataSource); 256 257 // 销毁pTargetGeometrys 258 OGRGeometryFactory::destroyGeometry(pTempGeometryUnion); 259 pTempGeometryUnion=NULL; 260 261 return true; 262 }