移动曲面法数字高程模型内插程序（C#桌面窗体）

运行环境：Visual Studio 2019/2022。别的不知道可不可以。

---

一、程序演示

DEM移动拟合插值算法成果演示

二、目的要求

1.算法要求

2.参考公式

移动拟合法插值的数学模型:

设计矩阵 M:  是一个 n×6 的矩阵，其中 n 是数据点的数量。每行代表一个数据点的坐标值映射到二次多项式项上。

权重矩阵 P:是一个对角矩阵，其对角线元素是每个数据点到目标点（x, y）的距离的倒数平方根的倒数。 

计算系数 ABCDEF: 使用加权最小二乘法，我们首先计算 ，然后求其逆 。最后通过计算出系数向量Z，其中 Z是坐标数据中的 Z值组成的向量。

具体算法请参考下面给出的代码，代码中有详细注释，应该都能看懂。

三、算法思路

这只是我的思路，不绝对。从最简单的代码开始写，成功运行后，一步步写的更复杂一些。

1.先写简易的版本，插值实现高程值Z的计算，在代码内部输入坐标数据。运算后会显示message弹窗。

2.删除上一步的坐标数据，固定地址读取excel文件，实现算法的运行。并且把excel中的内容在dataview控件中呈现出来。

3.将定死的待算坐标高程xy改为手动输入textBox中任意值。

4.删除固定地址，实现任意表格的导入，通过读取textbox中路径，任意打开excel表格用LinearRegression类中的函数计算里面的数据，得到ABCDEF系数。

5. 运算后textZ文本框会显示Z值

6.算法加入P（P为未知点距离已知点的距离平方的倒数）和M矩阵。（因为一开始的算法我写错了，后来检查出来后，又修改算法。这一步放到第一步中。）

7.完善代码，自定义输入xy坐标等功能

四、算法内容

这是我的源代码文件：

移动曲面法数字高程模型内插算法（C#桌面窗体版）资源-CSDN文库

1.安装程序包

mathnet.numerics包

具体看前辈文章：VisualStudio2022为C#安装mathnet.numerics包_mathnet.numerics下载安装-CSDN博客

好像还需要安装其他程序包，因为时间太久远，我记不清了，大家运行代码的时候，直接复制错误信息在CSDN搜索就能搜索出来需要安装什么包。

欢迎在评论区补充需要安装的包。

2.引用类库

3.具体代码

Form.cs

窗体布局

代码

using System;
using System.Collections.Generic;
using System.ComponentModel;
using System.Data;
using System.Drawing;
using System.Linq;
using System.Text;
using System.Threading.Tasks;
using System.Windows.Forms;
using System.Data.OleDb;
using System.Reflection;
using System.IO;
using MathNet.Numerics.LinearAlgebra;
using MathNet.Numerics.LinearAlgebra.Double;
using NPOI.SS.UserModel;
using NPOI.XSSF.UserModel; // 对于.xlsx文件  
using NPOI.HSSF.UserModel; // 对于.xls文件  
using static System.Windows.Forms.VisualStyles.VisualStyleElement;
using DEM插值算法;

namespace DEM插值算法
{
    public partial class Form1 : Form
    {
        private const int MinimumDataPoints = 6; // 假设最小数据点数量是6
                                                  // 将 data 提升为类级别的字段  
        public Form1()
        {
            InitializeComponent();
            comboBox1.Items.Add("移动曲面法");
        }

        //定义方法，读取Excel文件，并返回一个包含三元组（Tuple）的列表，每个三元组包含三个 double 类型的值
        private List<Tuple<double, double, double>> ReadExcelFile(string filePath)
        {
            List<Tuple<double, double, double>> data = new List<Tuple<double, double, double>>();

            using (FileStream stream = new FileStream(filePath, FileMode.Open, FileAccess.Read))
            {
                // 创建一个XSSFWorkbook实例，用于读取XLSX格式的Excel文件  
                IWorkbook workbook = new XSSFWorkbook(stream); 
                ISheet sheet = workbook.GetSheetAt(0); // 获取第一个工作表  

                if (sheet == null)
                {
                    throw new InvalidOperationException("Excel文件中没有找到工作表。");
                }

                IRow headerRow = sheet.GetRow(0); // 获取标题行

                // 从第二行开始读取数据  
                for (int rowNum = 1; rowNum <= sheet.LastRowNum; rowNum++)
                {
                    IRow row = sheet.GetRow(rowNum);   // 获取当前行  
                    
                    // 如果行不为空，则读取数据  
                    if (row != null) 
                    {
                        double x = 0.0, y = 0.0, z = 0.0;

                        // 获取x、y、z所在列的单元格  
                        ICell cellX = row.GetCell(0); // 假设x在第一列  
                        ICell cellY = row.GetCell(1); // 假设y在第二列  
                        ICell cellZ = row.GetCell(2); // 假设z在第三列  

                        // 如果单元格不为空，则将其转换为double类型并赋值给相应的变量
                        if (cellX != null) x = ConvertCellToDouble(cellX);
                        if (cellY != null) y = ConvertCellToDouble(cellY);
                        if (cellZ != null) z = ConvertCellToDouble(cellZ);

                        // 将读取到的x、y、z值添加到元组列表中
                        data.Add(Tuple.Create(x, y, z));
                    }
                }
            }
            return data;
        }

        //将Apache POI库中的ICell对象转换为double类型，否则无法正常运算
        private double ConvertCellToDouble(ICell cell)
        {
            switch (cell.CellType)
            { 
                case CellType.Numeric:// 如果单元格是数值类型  
                    return cell.NumericCellValue;// 直接返回单元格的数值  
                case CellType.String: // 如果单元格是字符串类型  
                    double result;  // 定义一个double类型的变量，用于存储转换结果  
                    if (double.TryParse(cell.StringCellValue, out result))// 尝试将字符串转换为double类型  
                    {
                        return result;
                    }
                    throw new InvalidOperationException("无法将字符串单元格转换为double类型。");
                default:
                    throw new InvalidOperationException("不支持的单元格类型。");
            }
        }

        private void CalculateCoefficientsButton_Click(object sender, EventArgs e)
        {
            string filePath = textBoxFilePath.Text; // 从textBoxFilePath中获取文件路径  

            // 确保文件路径不是空的  
            if (string.IsNullOrWhiteSpace(filePath))
            {
                MessageBox.Show("请输入有效的文件路径。");
                return;
            }

            try
            {
                // 1. 从Excel文件中读取数据  
                List<Tuple<double, double, double>> data = ReadExcelFile(filePath);

                // 确保读取到了足够的数据点来进行计算  
                if (data.Count < MinimumDataPoints)
                {
                    throw new InvalidOperationException("数据点数量不足以进行计算。");
                }
               if (comboBox1.SelectedIndex == 0)//根据下拉框执行移动拟合算法
               
               {
                    // 2. 调用 LinearRegression 类中的 ComputeABCDE 方法求解  
                    double x1 = double.Parse(textX.Text);
                    double y1 = double.Parse(textY.Text);
                    DenseMatrix coefficients = LinearRegression.ComputeABCDEF(data, x1, y1);

                    // 3. 使用计算得到的系数来计算新的 z 值  
                    // 从两个文本框中获取 x 和 y 的值  
                    if (double.TryParse(textX.Text, out double x) && double.TryParse(textY.Text, out double y))
                    {
                        // 使用 x 和 y 的值计算 z  
                        double z = LinearRegression.CalculateZ(coefficients, x, y);
                      
                        // 4. 输出结果  
                        MessageBox.Show($"对于 x={x} 和 y={y}, z={z}");
                    textZ.Text = $"{z}";
                    }
                    else
                    {
                    // 如果输入的不是有效的数字，显示错误消息  
                        MessageBox.Show("请输入有效的 x 和 y 值。");
                    }
               }
            }
            catch (Exception ex)  // 显示异常信息  
            {
               MessageBox.Show("发生错误：" + ex.Message);
            }
        }

        //打开文件对话框，选择excel文件，读取该文件的内容后，将数据展示在DataGridView控件中
        private void LoadExcelData_Click(object sender, EventArgs e)
        {
            using (OpenFileDialog openFileDialog = new OpenFileDialog())
            {

                openFileDialog.Title = "已知数据点坐标";
                openFileDialog.Filter = "高程点文件|*csv;*.xls;*.xlsx";
                if (openFileDialog.ShowDialog() == DialogResult.OK)
                {
                    string filePath = openFileDialog.FileName;
                    string connectionString = GetConnectionString(filePath);
                    textBoxFilePath.Text = openFileDialog.FileName;

                    // 使用OleDbConnection对象来建立与Excel文件的连接
                    using (OleDbConnection connection = new OleDbConnection(connectionString))  
                    {
                        try
                        {
                            connection.Open();

                            // 创建一个OleDbDataAdapter对象，用于执行SQL查询并填充数据集  
                            OleDbDataAdapter adapter = new OleDbDataAdapter("SELECT * FROM [Sheet1$]", connection);
                            DataSet dataSet = new DataSet();      // 创建一个新的数据集  
                            adapter.Fill(dataSet);   // 使用适配器填充数据集  
                            dataGridView1.DataSource = dataSet.Tables[0];// 将数据集的第一个表绑定到dataGridView1控件上 
                        }
                        catch (Exception ex)
                        {
                            MessageBox.Show("Error loading Excel data: " + ex.Message);
                        }
                    }
                }


            }
        }
     
        private string GetConnectionString(string filePath)
        {
            string connectionString = "";

            if (filePath.EndsWith(".xls"))
            {
                // Excel 97-2003 (.xls)  
                connectionString = "Provider=Microsoft.Jet.OLEDB.4.0;Data Source=" + filePath + ";Extended Properties='Excel 8.0;HDR=YES;';";
            }
            else if (filePath.EndsWith(".xlsx"))
            {
                // Excel 2007 and above (.xlsx)  
                connectionString = "Provider=Microsoft.ACE.OLEDB.12.0;Data Source=" + filePath + ";Extended Properties='Excel 12.0 Xml;HDR=YES;';";
            }

            return connectionString;
        }
    }
}



      
       

LinearRegression.cs

using MathNet.Numerics.LinearAlgebra;
using MathNet.Numerics.LinearAlgebra.Double;
using System;
using System.Collections.Generic;
using System.Linq;
using MathNet.Numerics;
using MathNet.Numerics.LinearAlgebra.Factorization;

namespace DEM插值算法
{
    internal class LinearRegression
    {
        //创建的对角权重矩阵P
        public static DenseMatrix ComputeWeightsMatrix(List<Tuple<double, double, double>> data, double x, double y)
        {
            var weights = data.Select(p => 1.0 / (Math.Pow(p.Item1 - x, 2) + Math.Pow(p.Item2 - y, 2))).ToList();

            //静态方法，用于计算线性回归的系数矩阵。它接收一个包含三元组的列表（每个三元组包含三个double类型的值）
            //以及两个double类型的参数x和y。这个方法返回一个DenseMatrix对象，表示系数矩阵。
            var P = DenseMatrix.OfDiagonalVector(DenseVector.OfEnumerable(weights));
            return P;
        }

        //创建M矩阵
        public static DenseMatrix ComputeMMatrix(List<Tuple<double, double, double>> data)
        {
            var M = DenseMatrix.OfArray(new double[data.Count, 6]);
            for (int i = 0; i < data.Count; i++)
            {
                var x = data[i].Item1;
                var y = data[i].Item2;
                M[i, 0] = x * x;      // X^2
                M[i, 1] = x * y;      // XY
                M[i, 2] = y * y;      // Y^2
                M[i, 3] = x;          // X
                M[i, 4] = y;          // Y
                M[i, 5] = 1.0;        // Constant
            }
            return M;
        }
        //计算ABCDEF系数
        public static DenseMatrix ComputeABCDEF(List<Tuple<double, double, double>> data, double x, double y)
        { 
            var P = ComputeWeightsMatrix(data, x, y);//调用ComputeWeightsMatrix方法计算权重矩阵P。
            var M = ComputeMMatrix(data);//调用ComputeMMatrix方法计算M矩阵。
            var MT = M.Transpose();//计算M矩阵的转置矩阵MT。

            var MP = MT * P * M;//计算MP矩阵，即MT、P和M的乘积。
            var MPInv = MP.Inverse();//计算MP矩阵的逆矩阵MPInv。

            if (MPInv == null)
            {
                throw new InvalidOperationException("Failed to compute inverse of matrix MP.");
            }

            var zValues = data.Select(p => p.Item3).ToArray();
            var MTP = MT * P;//计算MTP矩阵，即MT和P的乘积
            var ABCDEF = MPInv * MTP * DenseVector.OfArray(zValues);//计算线性回归的系数向量ABCDEF。
            // 创建一个 1 行 6 列的 DenseMatrix  
            var coefficientsMatrix = DenseMatrix.OfRowArrays(new[] { ABCDEF.ToArray() }); //将系数向量ABCDEF转换为一个1行6列的DenseMatrix对象。

            return coefficientsMatrix;
          
        }

        //计算待定的Z值
        public static double CalculateZ( DenseMatrix coefficients,double x, double y)
        {
            if (coefficients.RowCount != 1 || coefficients.ColumnCount != 6)
            {
                throw new ArgumentException("Coefficients matrix must be a 1x6 matrix.");
            }

            double A = coefficients[0, 0];
            double B = coefficients[0, 1];
            double C = coefficients[0, 2];
            double D = coefficients[0, 3];
            double E = coefficients[0, 4];
            double F = coefficients[0, 5];

            return A * x * x + B * x * y + C * y * y + D * x + E * y + F;
        }
    }
}

五、算法心得

应当由简入繁，在最简单代码的基础上一步步丰富内容，增加功能，变得复杂高级。

每一步都要运行出来后再进行下一步，不然全写完后再运行，代码太多根本不知道哪出问题了。