对于SuperMap的SFC（SuperMap Foundation Class）产品剖面分析很简单，只需要在栅格数据集（DEM或者GRID）上绘制一条线就可以直接取到这个剖面线上点的高程并绘制一条剖面线出来。对于SuperMap的UGC（Universal GIS Class）产品，并没有提供相应的接口，所以就自己模拟写了一个剖面分析的实现。
它山之石攻玉——确定实现的思路
1、数据准备
对照SuperMap Deskpro,实现剖面分析，需要一个栅格数据集和一条剖面线。通过SuperMap Objects.NET的Action.CreatePolyline或者SuperMap.UI.Action.CreateLine,可以得到一条剖面线；通过SuperMap.Data.Datasource.Datasets["name"/iIndex]获取一个栅格数据集，或者SuperMap.UI.MapControl.Map.Layers["name"/iIndex],获取栅格数据集对应的图层，然后Layer.Dataset获取栅格数据集。代码如下：
private void tsmiCreateProfile_Click(object sender, EventArgs e)
{
try
{
Layers Lyrs = null;
Lyrs = mapControl.Map.Layers;
if (Lyrs.Count > 0)
{
//判断是否有选中的线对象
if (IsNotEmptySelection())
{
SectionMap frm = new SectionMap();
frm.Show();
}
else
{
foreach (Layer Lyr in Lyrs)
{
Dataset dt = Lyr.Dataset;
if (dt.Type == DatasetType.Grid)
{
mapControl.TrackMode = SuperMap.UI.TrackMode.Track;
mapControl.TrackingStyle.LineColor = Color.Red;
mapControl.TrackingStyle.LineWidth = 0.2;
//mapControl.Action = SuperMap.UI.Action.CreateLine;
mapControl.Action = SuperMap.UI.Action.CreatePolyline;
bCreateProfile = true;
}
dt.Close();
}
}
}
}
catch (Exception ex)
{
MessageBox.Show("剖面分析出错:" + ex.Message);
throw;
}
}
/// <summary>
/// 返回地图中是否有选中的线对象，且选中对象的个数为1
/// </summary>
/// <returns>True Or False</returns>
private Boolean IsNotEmptySelection()
{
Boolean bSelectOneLine = false;
try
{
Selection[] sels = mapControl.Map.FindSelection(true);
foreach (Selection sel in sels)
{
Dataset selDt = sel.Dataset;
if (selDt.Type == DatasetType.Line && sel.Count == 1)
{
bSelectOneLine = true;
Recordset rec = sel.ToRecordset();
ProfileLine = rec.GetGeometry() as GeoLine;
rec.Close();
rec.Dispose();
sel.Dispose();
break;
}
else
{
sel.Dispose();
continue;
}
}
return bSelectOneLine;
}
catch (Exception ex)
{
MessageBox.Show("剖面分析出错:" + ex.Message);
return bSelectOneLine;
}
}
private void mapControl_Tracked(object sender, TrackedEventArgs e)
{
Geometry geo = e.Geometry;
try
{
if (geo.IsEmpty || geo.Type != GeometryType.GeoLine )
{
return;
}
else
{
Layers Lyrs = null;
Lyrs = mapControl.Map.Layers;
if (Lyrs.Count > 0)
{
foreach (Layer Lyr in Lyrs)
{
Dataset dt = Lyr.Dataset;
if (dt.Type == DatasetType.Grid)
{
SectionMap frm = new SectionMap(dt as DatasetGrid, geo as GeoLine);
frm.Show();
}
dt.Close();
}
}
}
workspaceControl.Refresh();
}
catch (Exception ex)
{
MessageBox.Show("出错:" + ex.Message);
}
}
2、算法准备
在《空间分析建模与原理》的183页介绍了剖面分析的相关算法，这里不再赘言。通过该算法需要得到剖面线和栅格格网数据集的交点并得到每个交点上的高程值或者像素值。
a、得到采样点
假定剖面线只是一条简单的直线段，可以计算出选段的起点终点间XY方向的增量dX和dY，根据增量的情况分以下四类情况计算采样点：
dX=0,此时剖面线和栅格格网的纵轴方向一致，采样点坐标为(a,yi),i = 1,2...int(dy/dresolutiony)(dresolutionY是栅格格网Y轴的分辨率)； dY=0,此时剖面线和栅格格网的横轴方向一致，采样点坐标为（xj，b）,j=1,2...int(dx/dresolutionx)(dresolutionX是栅格格网X轴的分辨率)； Abs(dy/dx)<=1,dx!=0,采样点的y坐标为：yk = (dy/dx)(xk-x1)+y1,k=1,2...int(dx/dresolutionx),采样点坐标为（xk,yk）; Abs(dy/dx)>1,dx!=0,采样点的X坐标为：xs = (dy/dx)(ys-y1)+x1,k=1,2...int(dy/dresolutiony)，采样点坐标为（xs,ys）;
b、得到采样点的栅格值
把采样点的坐标转换成栅格格网的行列号，通过SuperMap.Data.DatasetGrid.XYToGrid； 通过行列号获取采样点的栅格值，SuperMap.Data.DatasetGrid.GetValue; 代码如下：public static System.Data.DataTable CreateResamplePoint2DEx(DatasetGrid dtg, Point2D ptStart, Point2D ptEnd, double dResolutionX, double dResolutionY, double dIncrementDistance) { //剖面线采样点 System.Data.DataTable myList = new System.Data.DataTable(); try { if (ptStart.IsEmpty || ptEnd.IsEmpty) { return myList; } else { //生成Tabable表 myList.Columns.Add("ID"); myList.Columns.Add("X"); myList.Columns.Add("Y"); myList.Columns.Add("Distance"); myList.Columns.Add("Z"); //起点终点间的增量 double dIncrementX = Math.Abs(ptStart.X - ptEnd.X); double dIncrementY = Math.Abs(ptStart.Y - ptEnd.Y); //线段上采样点的个数 int iPtsCount = 0; #region X方向增量为0 //第一種情况：X方向增量为0 if (dIncrementX == 0) { iPtsCount = Convert.ToInt32(dIncrementY / dResolutionY); for (int i = 0; i <= iPtsCount; i++) { //采样点的平面坐标 Point2D pt2D = new Point2D(); pt2D.X = ptStart.X; pt2D.Y = ptStart.Y + dResolutionY * i; //采样点到起点的距离 double dDistance = Math.Sqrt(Math.Pow(pt2D.X - ptStart.X, 2) + Math.Pow(pt2D.Y - ptStart.Y, 2)); //采样点的高程,判断坐标是否不再栅格数据集的地理范围之内，如果在就计算栅格值，否则跳过不计算 double dElevation = 0; if (BeyondBounds(dtg, pt2D)) { //采样点所对应的行列号 Point ptRowCol = dtg.XYToGrid(pt2D); dElevation = dtg.GetValue(ptRowCol.X, ptRowCol.Y); } else { dElevation = 0; } //添加采样点 DataRow dr = myList.NewRow(); dr[0] = i; dr[1] = pt2D.X; dr[2] = pt2D.Y; dr[3] = dDistance + dIncrementDistance; dr[4] = dElevation; myList.Rows.Add(dr); } } //第二種情况：Y方向增量为0 else if (dIncrementY == 0) { iPtsCount = Convert.ToInt32(dIncrementX / dResolutionX); for (int i = 0; i < iPtsCount; i++) { //采样点的平面坐标 Point2D pt2D = new Point2D(); pt2D.X = ptStart.X + dResolutionX * i; pt2D.Y = ptStart.Y; //采样点到起点的距离 double dDistance = Math.Sqrt(Math.Pow(pt2D.X - ptStart.X, 2) + Math.Pow(pt2D.Y - ptStart.Y, 2)); //采样点的高程,判断坐标是否不再栅格数据集的地理范围之内，如果在就计算栅格值，否则跳过不计算 double dElevation = 0; if (BeyondBounds(dtg, pt2D)) { //采样点所对应的行列号 Point ptRowCol = dtg.XYToGrid(pt2D); dElevation = dtg.GetValue(ptRowCol.X, ptRowCol.Y); } else { dElevation = 0; } //添加采样点 DataRow dr = myList.NewRow(); dr[0] = i; dr[1] = pt2D.X; dr[2] = pt2D.Y; dr[3] = dDistance + dIncrementDistance; dr[4] = dElevation; myList.Rows.Add(dr); } } #endregion #region X方向增量不为0，线段有斜率的情况 double dSlope = dIncrementY / dIncrementX; //第三种種情况：X方向增量不为0，且斜率小于等于1 if (dIncrementX != 0 && dSlope <= 1) { iPtsCount = Convert.ToInt32(dIncrementX / dResolutionX); for (int i = 0; i < iPtsCount; i++) { //采样点的平面坐标 Point2D pt2D = new Point2D(); pt2D.X = ptStart.X + dResolutionX * i; pt2D.Y = dSlope * (pt2D.X - ptStart.X) + ptStart.Y; //采样点到起点的距离 double dDistance = Math.Sqrt(Math.Pow(pt2D.X - ptStart.X, 2) + Math.Pow(pt2D.Y - ptStart.Y, 2)); //采样点的高程,判断坐标是否不再栅格数据集的地理范围之内，如果在就计算栅格值，否则跳过不计算 double dElevation = 0; if (BeyondBounds(dtg, pt2D)) { //采样点所对应的行列号 Point ptRowCol = dtg.XYToGrid(pt2D); dElevation = dtg.GetValue(ptRowCol.X, ptRowCol.Y); } else { dElevation = 0; } //添加采样点 DataRow dr = myList.NewRow(); dr[0] = i; dr[1] = pt2D.X; dr[2] = pt2D.Y; dr[3] = dDistance + dIncrementDistance; dr[4] = dElevation; myList.Rows.Add(dr); } } //第四种種情况：X方向增量不为0，且斜率大于1 else if (dIncrementX != 0 && dSlope > 1) { iPtsCount = Convert.ToInt32(dIncrementY / dResolutionY); for (int i = 0; i < iPtsCount; i++) { //采样点的平面坐标 Point2D pt2D = new Point2D(); pt2D.Y = ptStart.Y + dResolutionY * i; pt2D.X = (pt2D.Y - ptStart.Y) / dSlope + ptStart.X; //采样点到起点的距离 double dDistance = Math.Sqrt(Math.Pow(pt2D.X - ptStart.X, 2) + Math.Pow(pt2D.Y - ptStart.Y, 2)); ; //采样点的高程,判断坐标是否不再栅格数据集的地理范围之内，如果在就计算栅格值，否则跳过不计算 double dElevation = 0; if (BeyondBounds(dtg, pt2D)) { //采样点所对应的行列号 Point ptRowCol = dtg.XYToGrid(pt2D); dElevation = dtg.GetValue(ptRowCol.X, ptRowCol.Y); } else { dElevation = 0; } //添加采样点 DataRow dr = myList.NewRow(); dr[0] = i; dr[1] = pt2D.X; dr[2] = pt2D.Y; dr[3] = dDistance + dIncrementDistance; dr[4] = dElevation; myList.Rows.Add(dr); } } #endregion //返回采样点距离列表类 return myList; } } catch (Exception ex) { MessageBox.Show("生成采样点数据出错:" + ex.Message); return myList; } }c、绘制剖面线
绘制剖面线有三种方法：借助MSChart控件、通过底层API、借助Excel。本文中使用第一種方法。
代码如下：
private void CreateChart() { if (dtgProfile == null || lineProfile.IsEmpty) { return; } else { //创建剖面线对象 createProfile = new CreateProfile(dtgProfile, lineProfile); //数据表对象 dataTable = createProfile.CreateResampleDataEx(dtgProfile, lineProfile); //计算高程的最大值和最小值 string dMaxZ = dataTable.Compute("Max(Z)","").ToString(); string dMinZ = dataTable.Compute("Min(Z)", "").ToString(); Series series = chartMap.Series[0]; series.ChartType = SeriesChartType.SplineArea; series.BorderWidth = 2; series.XAxisType = AxisType.Primary; series.XValueType = ChartValueType.Auto; series.XValueMember = "Date"; series.YAxisType = AxisType.Primary; series.YValueType = ChartValueType.Int32; series.YValuesPerPoint = 1; series.YValueMembers = "UnitsSold"; series.IsXValueIndexed = true; chartMap.Titles.Add("MaxZ = " + dMaxZ + " MinZ = " + dMinZ); foreach (DataRow row in dataTable.Rows) { double dDistance = Convert.ToDouble(row["Distance"]); int iDistance = Convert.ToInt32(Math.Ceiling(dDistance)); series.Points.AddXY(iDistance, row["Z"]); } chartMap.ChartAreas[0].Area3DStyle.Enable3D = true; //展示 3D } }
结果展示
不足说明：
由于能力所限，不免有所疏漏，剖面线显得不够光滑。希望各位网友能给予指正，谢谢。
