基于gdal 和 shp求解多点最短路径

大致思路是通过将多点进行排列组合，求解出来地理位置上的最短路径，在通过拓扑构网，将shp文件构建成图，通过求解两点之间的最优解来求解最短路径
全排求解最短路径

class TDPermutations
{
public:
	std::vector<OGRPoint>  calPermutations(std::vector<OGRPoint> vecOGRPoint);
private:
	std::vector<OGRPoint> m_vecOGRPoint;
	std::vector<std::vector<OGRPoint>> m_indexandPointInfo;
	std::vector<std::pair<int,double>> m_pathLang;
};

std::vector<OGRPoint> TDPermutations::calPermutations(std::vector<OGRPoint> vecOGRPoint)
{
	m_vecOGRPoint = vecOGRPoint;
	do {
		std::vector<OGRPoint> temp = m_vecOGRPoint;
		m_indexandPointInfo.push_back(temp);
	} while (std::next_permutation(m_vecOGRPoint.begin(), m_vecOGRPoint.end(), [](OGRPoint & it1, OGRPoint &it2)->bool{
		return it1.getX() == it2.getX() ? it1.getY() < it2.getY() : it1.getX() <it2.getX();
	}));
	int num = 0;
	std::for_each(m_indexandPointInfo.begin(), m_indexandPointInfo.end(), [&](std::vector<OGRPoint> &temp){
		double linelang = 0.0;
		for (auto index = 0; index < temp.size() - 1; index++)
		{
			linelang += temp[index].Distance(&temp[index + 1]);
		}
		m_pathLang.push_back({ num,linelang });
		num++;
	});
	auto it = std::min_element(m_pathLang.begin(), m_pathLang.end(), [](std::pair<int,double> & it1, std::pair<int, double> & it2)->bool{
		if (it1.second < it2.second)
		{
			return true;
		}
		return false;
	});
	return m_indexandPointInfo[it->first];
}

shp转化成GNMGraph对象求解两点之间的最短路径

class TDOGRPoint :public OGRPoint
{
public:
	bool  operator < (const TDOGRPoint& str) const
	{ 
		return getX() == str.getX() ? getY() < str.getY() : getX() < str.getX();
	 
	}
	const bool  operator > (const TDOGRPoint& str) const
	{
		return getX() == str.getX() ? getY() < str.getY() : getX() < str.getX();
	}
	bool operator == (const TDOGRPoint& str) const
	{
		return getX() == str.getX() && getY() == str.getY();
	}
};
class TDGdalRoadgraph :public TDconversionBase{
public:
	void process(LayerPathInfo, LayerPathInfo);
	OGRPoint beginPoint;
	OGRPoint endPoint;
private:
	bool readOriginalInfo(LayerPathInfo originalPathInfo, LayerPathInfo targetpatahinfo);
	void makeGraph();
	void tranformToPoint(TranForm_enum enum_obj, OGRFeature *poFeature) { }
private:
	std::vector<OGRPoint> m_params;
	std::map<TDOGRPoint,int> m_mapPointindex;
	std::map<int, TDOGRPoint> m_mapindexPoint;
	std::vector<OGRLineString> m_tempOGRLine;
	GNMGraph m_nmGraph;
	int beginFID = { 0 };
	int endFID = { 0 };
	double beginlangPath = { 0.0 };
	double endlangPath = { 0.0 };
};

void TDGdalRoadgraph::process(LayerPathInfo sourpathinfo, LayerPathInfo targetpatahinfo)
{
	readOriginalInfo(sourpathinfo, targetpatahinfo);
}
void TDGdalRoadgraph::makeGraph()
{
	auto fun2 = [&](double &ogrlang, OGRPoint & it, OGRPoint & it1,int innum,int & autnum)
	{
		if (OGR_G_Distance(&it, &it1) < ogrlang)
		{
			ogrlang = OGR_G_Distance(&it, &it1);
			autnum = innum;
		}
	};
	auto fun1 = [&](OGRPoint & it,int num) {
		fun2(beginlangPath, beginPoint, it, num, beginFID);
		fun2(endlangPath, endPoint, it, num, endFID);
	};
	if (!m_orgLayer)
		return;
	OGRFeature * feature = NULL;
	m_orgLayer->ResetReading();
	int FDID = 0;
	int lineFDID = 0;
	while ((feature = m_orgLayer->GetNextFeature()) != NULL)
	{
		OGRMultiLineString  mpl;
		if (feature->GetGeometryRef()->getGeometryType() == wkbMultiLineString)
			mpl = (*feature->GetGeometryRef()->toMultiLineString());
		else if (feature->GetGeometryRef()->getGeometryType() == wkbLineString)
			mpl.addGeometry(feature->GetGeometryRef()->toLineString());
		auto geomNum = mpl.getNumGeometries();
		for (int i = 0; i < geomNum; i++)
		{
			OGRPoint pt; TDOGRPoint pts;
			auto lineString = mpl.getGeometryRef(i)->toLineString();
			auto numpoint = lineString->getNumPoints();
			for (int j = 0; j < numpoint; j++)
			{
				lineString->getPoint(j, &pt);
				pts.setX(pt.getX());
				pts.setY(pt.getY());
				if (m_mapPointindex.find(pts) == m_mapPointindex.end())
				{
					m_mapPointindex[pts] = FDID;
					m_mapindexPoint[FDID] = pts;
					if (FDID == 0)
					{
						beginlangPath = OGR_G_Distance(&beginPoint, &pt);
						endlangPath = OGR_G_Distance(&endPoint, &pt);
					}
					fun1(pt, FDID);
					m_nmGraph.AddVertex(FDID);
					FDID++;
				}
			}
			for (int j = 0; j < numpoint -1; j++)
			{
				OGRLineString ogrlinestring;
				OGRPoint *pt1 = new OGRPoint(), *pt2 = new OGRPoint();
				lineString->getPoint(j, pt1);
				ogrlinestring.addPoint(pt1);
				lineString->getPoint(j+1, pt2);
				ogrlinestring.addPoint(pt2);
				m_tempOGRLine.push_back(ogrlinestring);
				lineFDID++;
			}
		}
	}//将节点存储到图
	//遍历数据,将道路长度和
	lineFDID = 0;
	std::for_each(m_tempOGRLine.begin(), m_tempOGRLine.end(), [&](OGRLineString & ogrlinestring) {
		OGRPoint pt,pt1;
		ogrlinestring.getPoint(0, &pt);
		TDOGRPoint pts;
		pts.setX(pt.getX());
		pts.setY(pt.getY());
		auto beginindex = m_mapPointindex[pts];
		ogrlinestring.getPoint(1, &pt1);
		pts.setX(pt1.getX());
		pts.setY(pt1.getY());
		auto endindex = m_mapPointindex[pts];
		m_nmGraph.AddEdge(lineFDID, beginindex, endindex, true, OGR_G_Distance(&pt, &pt1), OGR_G_Distance(&pt, &pt1));
		lineFDID++;
	});
}
bool TDGdalRoadgraph::readOriginalInfo(LayerPathInfo originalPathInfo, LayerPathInfo targetpatahinfo)
{
	(void)openOriginalInfo(originalPathInfo);
	IF_NULL_RETRUN_FALSE(m_orgLayer)
		OGRFeature *poFeature;
	m_orgLayer->ResetReading();
	OGRFeatureDefn * ogrfeaturedefn = m_orgLayer->GetLayerDefn();
	pmoSourceSRS = m_orgLayer->GetSpatialRef();
	makeGraph();
	auto indexinfo = m_nmGraph.DijkstraShortestPath(beginFID,endFID);
	std::for_each(indexinfo.begin(), indexinfo.end(), [&](EDGEVERTEXPAIR & indexit) {
		
		if (indexit.first != -1 && indexit.second != -1)
		{
			OGRLineString  ogrlinestring;
			ogrlinestring = m_tempOGRLine[indexit.second];
			OGRFeatureDefn * ogrfeaturedefn = new OGRFeatureDefn();
			OGRFeature * ogrfeature = new OGRFeature(ogrfeaturedefn);
			ogrfeature->SetGeometry(&ogrlinestring);
			m_cachevector.push_back(ogrfeature);
		}
	});
	wiriteInfo(targetpatahinfo, ogrfeaturedefn);
	m_cachevector.clear();
	if (NULL != m_pwriteDataset)
	{
		GDALClose(m_pwriteDataset);
	}
	if (NULL != m_preadDataset)
	{
		GDALClose(m_preadDataset);
	}
	return true;
}