当前位置：首页 > news >正文

【QGIS二次开发】空间分析-10

news 来源：原创 2025/6/4 5:20:08

系列目录：

【QGIS二次开发】地图显示与交互-01_qgis二次开发加载地图案例-CSDN博客

【QGIS二次开发】地图显示与交互-02_setlayerlabeling-CSDN博客

【QGIS二次开发】地图符号与色表-03-CSDN博客

【QGIS二次开发】地图编辑-04-CSDN博客

【QGIS二次开发】地图编辑-05-CSDN博客

【QGIS二次开发】地图编辑-06-CSDN博客

【QGIS二次开发】地图编辑-07-CSDN博客

【QGIS二次开发】地图编辑-08-CSDN博客

【QGIS二次开发】地图编辑-09-CSDN博客

5.1 矢量缓冲区分析

功能点描述：

矢量缓冲区分析：分别实现点、线路和区的缓冲区，并输出结果区文件。

运行截图：

操作前：

图 89 示例数据

操作后：

图 90 缓冲区结果

主要代码如下：

该代码是一个执行缓冲区分析的函数performBufferAnalysis。首先，它通过GDAL库注册所有驱动器，打开输入文件，获取输入图层，并创建输出文件和输出图层。接着，它使用OGRGeometryFactory创建一个用于缓冲区分析的几何对象(wkbPolygon)。然后，通过遍历输入图层的要素，获取每个要素的几何对象，并对几何对象执行缓冲区分析。缓冲距离和段数通过参数传递，生成缓冲区分析后的几何对象。

随后，创建输出要素，将缓冲区分析后的几何对象设置为输出要素的几何对象，并将输出要素写入输出图层。最后，释放内存，关闭输入文件和输出文件，弹出提示框告知用户缓冲区分析已完成。然后，通过QgsVectorLayer类将输出文件作为图层添加到QgsProject中，并刷新地图画布，以在地图中显示缓冲区分析的结果。

void ZoneAnalysis::performBufferAnalysis(const QString& inputFilePath, const QString& outputFilePath, double bufferDistance) {    // 注册所有的驱动器    GDALAllRegister();    // 打开输入文件    GDALDataset* inputDataset = (GDALDataset*)GDALOpenEx(inputFilePath.toStdString().c_str(), GDAL_OF_VECTOR, nullptr, nullptr, nullptr);    if (inputDataset == nullptr) {    QMessageBox::critical(this, "错误", "无法打开输入文件！");    return;    }    // 获取输入图层    OGRLayer* inputLayer = inputDataset->GetLayer(0);    // 创建输出文件    GDALDriver* driver = GetGDALDriverManager()->GetDriverByName("ESRI Shapefile");    GDALDataset* outputDataset = driver->Create(outputFilePath.toStdString().c_str(), 0, 0, 0, GDT_Unknown, nullptr);    if (outputDataset == nullptr) {    QMessageBox::critical(this, "错误", "无法创建输出文件！");    GDALClose(inputDataset);    return;    }    // 创建输出图层    OGRLayer* outputLayer = outputDataset->CreateLayer("buffer", nullptr, wkbPolygon, nullptr);    if (outputLayer == nullptr) {    QMessageBox::critical(this, "错误", "无法创建输出图层！");    GDALClose(outputDataset);    GDALClose(inputDataset);    return;    }    // 创建缓冲区分析对象    OGRGeometryFactory::createGeometry(wkbPolygon);    // 遍历输入图层的要素    inputLayer->ResetReading();    OGREnvelope envelope;    while (OGRFeature* feature = inputLayer->GetNextFeature()) {    OGRGeometry* geometry = feature->GetGeometryRef();    if (geometry != nullptr) {    // 执行缓冲区分析    OGRGeometry* bufferGeometry = geometry->Buffer(bufferDistance, 30);    // 创建输出要素    OGRFeature* outputFeature = OGRFeature::CreateFeature(outputLayer->GetLayerDefn());    outputFeature->SetGeometry(bufferGeometry);    // 将输出要素写入输出图层    outputLayer->CreateFeature(outputFeature);    // 释放内存    OGRFeature::DestroyFeature(outputFeature);    OGRGeometryFactory::destroyGeometry(bufferGeometry);    }    // 释放内存    OGRFeature::DestroyFeature(feature);    }    // 释放资源    GDALClose(outputDataset);    GDALClose(inputDataset);    // 缓冲区分析完成后，关闭缓冲区框    this->close();    // 弹出提示框告知用户    QMessageBox::information(this, "提示", "缓冲区分析已完成！");    // 从输出路径打开shp    QFileInfo fi(outputFilePath);    if (!fi.exists()) { return; }    QString layerBaseName = fi.baseName(); // 图层名称    QgsVectorLayer* vecLayer = new QgsVectorLayer(outputFilePath, layerBaseName);    if (!vecLayer) { return; }    if (!vecLayer->isValid())    {    QMessageBox::information(0, "", "layer is invalid");    return;    }    QgsProject::instance()->addMapLayer(vecLayer);    m_mapCanvas->refresh();    //m_mapCanvas->setExtent(vecLayer->extent());    //layers.append(vecLayer);    //m_mapCanvas->setLayers(layers);    //m_mapCanvas->refresh();    
}

5.2 矢量裁剪分析

功能点描述：

矢量裁剪分析：分别用矩形、圆和任意多边形对矢量图层进行裁剪，并将裁剪结果输出到矢量文件。

运行截图：

操作前：

图 91 裁剪数据

操作后：

图 92 裁剪效果

主要代码：

这段代码是一个名为ZoneCut的类的成员函数，用于对两个矢量文件进行裁剪，并将裁剪后的结果输出到一个新的矢量文件中。首先，它初始化GDAL库。然后，通过打开两个输入文件，创建输出文件，以及获取相应的图层等操作，准备进行矢量文件的裁剪操作。如果打开输入文件或创建输出文件失败，则函数会提前返回。

在裁剪操作的过程中，它遍历第一个输入文件的每一个特征，获取其几何体，并为每个特征创建一个裁剪几何体的克隆。接着，将第二个输入文件的空间过滤器设置为裁剪几何体，遍历第二个输入文件的每一个特征，如果其几何体与裁剪几何体相交，则计算它们的交集，并将交集作为输出图层的一个特征。

最后，完成裁剪操作后，关闭输入文件和输出文件，清空输入和输出文件路径，同时关闭与裁剪相关的对话框。该函数主要实现了两个矢量文件的裁剪，将裁剪后的结果保存到一个新的矢量文件中，为用户提供了对地理数据进行空间裁剪的功能。

void ZoneCut::performZoneCut(const QString& inputFilePath1, const QString& inputFilePath2, const QString& outputFilePath)    
{    GDALAllRegister();    // 打开第一个输入文件    QString inputFilePath1Utf8 = inputFilePath1.toUtf8();    GDALDataset* inputDataset1 = (GDALDataset*)GDALOpenEx(inputFilePath1.toStdString().c_str(), GDAL_OF_VECTOR, nullptr, nullptr, nullptr);    if (inputDataset1 == nullptr)    {    // 处理打开文件失败的情况    return;    }    // 打开第二个输入文件    GDALDataset* inputDataset2 = (GDALDataset*)GDALOpenEx(inputFilePath2.toStdString().c_str(), GDAL_OF_VECTOR, nullptr, nullptr, nullptr);    if (inputDataset2 == nullptr)    {    // 处理打开文件失败的情况    GDALClose(inputDataset1);    return;    }    // 创建输出文件    GDALDriver* outputDriver = GetGDALDriverManager()->GetDriverByName("ESRI Shapefile");    GDALDataset* outputDataset = outputDriver->Create(outputFilePath.toStdString().c_str(), 0, 0, 0, GDT_Unknown, nullptr);    if (outputDataset == nullptr)    {    // 处理创建文件失败的情况    GDALClose(inputDataset1);    GDALClose(inputDataset2);    return;    }    // 获取输入文件的第一个图层    OGRLayer* inputLayer1 = inputDataset1->GetLayer(0);    // 获取输入文件的第二个图层    OGRLayer* inputLayer2 = inputDataset2->GetLayer(0);    // 创建输出图层    OGRLayer* outputLayer = outputDataset->CreateLayer("result", inputLayer1->GetSpatialRef()->Clone(), inputLayer1->GetLayerDefn()->GetGeomType(), nullptr);    // 获取第一个输入图层的几何体    OGRFeature* inputFeature1 = nullptr;//初始化输入文件1的特征指针。    while ((inputFeature1 = inputLayer1->GetNextFeature()) != nullptr)    {    OGRGeometry* inputGeometry1 = inputFeature1->GetGeometryRef();    OGRGeometry* cutGeometry = inputGeometry1->clone();//克隆输入文件1的特征的几何体    // 设置第二个输入图层的空间过滤器为被裁剪几何体    inputLayer2->SetSpatialFilter(cutGeometry);    // 获取第二个输入图层的几何体    OGRFeature* inputFeature2 = nullptr;    while ((inputFeature2 = inputLayer2->GetNextFeature()) != nullptr)    {    OGRGeometry* inputGeometry2 = inputFeature2->GetGeometryRef();// 获取输入文件2的特征的几何体。    if (inputGeometry2->Intersects(cutGeometry))//如果输入文件2的特征的几何体与裁剪几何体相交。    {    OGRFeature* outputFeature = OGRFeature::CreateFeature(outputLayer->GetLayerDefn());    outputFeature->SetGeometry(inputGeometry2->Intersection(cutGeometry));//设置输出图层的特征的几何体为输入文件2的特征的几何体与裁剪几何体的交集。    outputLayer->CreateFeature(outputFeature);    OGRFeature::DestroyFeature(outputFeature);    }    OGRFeature::DestroyFeature(inputFeature2);    }    OGRFeature::DestroyFeature(inputFeature1);    OGRGeometryFactory::destroyGeometry(cutGeometry);    }    // 关闭数据集    GDALClose(inputDataset1);    GDALClose(inputDataset2);    GDALClose(outputDataset);    // 关闭裁剪对话框    accept();    ui.InputlineEdit1->clear();    ui.InputlineEdit2->clear();    ui.OutputlineEdit->clear();    
}

5.3 矢量转栅格

功能点描述：将矢量图转化为栅格图层显示挂树，并输出tiff文件。

操作及截图：

当点击输入图层按钮时，弹出文件选择对话框，用户可以选择输入的Shapefile文件。选择后，将文件路径显示在界面上的InputEdit文本框中。然后使用GDAL库打开选定的Shapefile文件，并获取第一个图层的字段名称，将这些字段名称添加到界面上的下拉框(comboBox)中供用户选择。

图 93 获取输入图层

选择属性的combobox是获取输入图层的属性字段作为矢量转换字段，之后点击输出输出图层按钮弹出文件保存对话框，用户可以选择输出的栅格文件（TIFF文件）的保存路径。选择后，将文件路径显示在界面上的OutputEdit文本框中。

图 94 选定属性字段，确定输出路径

最后点击确定按钮时，开始进行转换操作。首先，获取输入和输出文件的路径以及用户选择的属性名称，注册GDAL驱动并打开输入Shapefile文件，获取输入Shapefile文件的第一个图层，并获取图层的空间范围和投影信息并创建一个栅格驱动，并创建一个指定大小的栅格数据集，用于保存转换后的栅格数据。根据图层的空间范围和大小，设置栅格数据集的地理转换参数（地理坐标系、像元大小等）。如果图层具有投影信息，将投影信息设置到栅格数据集中，设置好后创建一个包含要转换图层的数组，并设置转换选项。使用GDAL的GDALRasterizeLayers函数将图层转换为栅格数据。最后打开转换后的栅格文件，并将其添加到QGIS项目中进行显示。

图 95 矢量转栅格成果

核心代码：

void shptranstiff::on_InputpushButton_clicked() {  QString filePath = QFileDialog::getOpenFileName(nullptr, "Select shapefile", "", "Shapefiles (*.shp)");  ui.InputEdit->setText(filePath);  GDALAllRegister();  // 注册GDAL驱动  // 打开输入数据集  GDALDataset* inputDataset = (GDALDataset*)GDALOpenEx(filePath.toUtf8().constData(), GDAL_OF_VECTOR, nullptr, nullptr, nullptr);  if (inputDataset == nullptr) {  QMessageBox::critical(this, "错误", "无法打开输入文件。");  return;  }  // 获取图层数量  int layerCount = inputDataset->GetLayerCount();  // 清空下拉框  ui.comboBox->clear();  // 获取选定图层索引  int selectedLayerIndex = 0; // 假设选择第一个图层  OGRLayer* selectedLayer = inputDataset->GetLayer(selectedLayerIndex);  // 获取字段数量  int fieldCount = selectedLayer->GetLayerDefn()->GetFieldCount();  // 用属性名称填充下拉框  for (int i = 0; i < fieldCount; i++) {  OGRFieldDefn* fieldDefn = selectedLayer->GetLayerDefn()->GetFieldDefn(i);  const char* attributeName = fieldDefn->GetNameRef();  ui.comboBox->addItem(QString::fromUtf8(attributeName));  }  // 关闭Shapefile  GDALClose(inputDataset);  
}  void shptranstiff::on_OutputpushButton_clicked() {  QString saveFilePath = QFileDialog::getSaveFileName(nullptr, "Save tiff file", "", "files (*.tiff)");  ui.OutputEdit->setText(saveFilePath);  
}  void shptranstiff::on_pushButton_clicked() {  QString inputPath = ui.InputEdit->text();  QString outputPath = ui.OutputEdit->text();  // 获取选定的属性名称  QString attributeName = ui.comboBox->currentText();  // 检查是否选择了属性名称  if (attributeName.isEmpty()) {  QMessageBox::critical(this, "错误", "请选择属性名称。");  return;  }  // 注册所有的驱动程序  GDALAllRegister();  // 打开矢量文件  GDALDataset* shpData = (GDALDataset*)GDALOpenEx(inputPath.toUtf8().constData(), GDAL_OF_VECTOR, nullptr, nullptr, nullptr);  if (shpData == nullptr) {  QMessageBox::critical(this, "Error", "Failed to open vector file.");  return;  }  // 获取第一个图层  OGRLayer* shpLayer = shpData->GetLayer(0);  OGREnvelope env;  shpLayer->GetExtent(&env);//获取图层的坐标范围到env指向的内存中  int m_nImageWidth = 1440;  int m_nImageHeight = 720;  OGRSpatialReference* pOgrSRS = shpLayer->GetSpatialRef();//数据投影信息  char* pPrj = NULL;  if (pOgrSRS == NULL)  {  QMessageBox::critical(this, "Error", "没有投影信息");  m_nImageHeight = (int)env.MaxX;  m_nImageWidth = (int)env.MaxY;  }  else  {  pOgrSRS->exportToWkt(&pPrj);  }  // 创建栅格驱动  GDALDriver* poDriver = GetGDALDriverManager()->GetDriverByName("GTiff");  GDALDataset* poNewDS = poDriver->Create(outputPath.toStdString().c_str(), m_nImageWidth, m_nImageHeight, 1, GDT_Float32, NULL);  double adfGeoTransform[6];  adfGeoTransform[0] = env.MinX;//左上角经度  adfGeoTransform[1] = (env.MaxX - env.MinX) / m_nImageWidth;//像元宽度(影像在宽度上的分辨率)  adfGeoTransform[2] = 0; //如果影像是指北的, 这个参数的值为0  adfGeoTransform[3] = env.MaxY;//左上角纬度  adfGeoTransform[4] = 0; //如果影像是指北的, 这个参数的值为0。  adfGeoTransform[5] = (env.MinY - env.MaxY) / m_nImageHeight; //像元高度(影像在高度上的分辨率)  GDALSetGeoTransform(poNewDS, adfGeoTransform);  if (pOgrSRS != NULL)  {  poNewDS->SetProjection(pPrj);  }  int* pnbandlist = new int[1];  pnbandlist[0] = 1;  OGRLayerH* player;  player = new OGRLayerH[1];  player[0] = (OGRLayerH)shpLayer;  char** papszOptions = NULL;  papszOptions = CSLSetNameValue(papszOptions, "CHUNKSIZE", "1");  papszOptions = CSLSetNameValue(papszOptions, "ATTRIBUTE", attributeName.toUtf8().constData());  CPLErr err = GDALRasterizeLayers((GDALDatasetH)poNewDS, 1, pnbandlist, 1, player,NULL, NULL, NULL, papszOptions, NULL, NULL);  GDALClose(shpData);  GDALClose(poNewDS);  GDALDestroyDriverManager();  delete[]player;  delete[]pnbandlist;  close();  QMessageBox::information(this, "成功", "栅格文件转换完成。");  // 打开转换后的栅格文件  QFileInfo fi(outputPath);  if (!fi.exists()) { return; }  QString layerBaseName = fi.baseName(); // 图层名称  QgsRasterLayer* rasterLayer = new QgsRasterLayer(outputPath, layerBaseName);  if (!rasterLayer->isValid())  {  QMessageBox::critical(this, "Error", "Failed to load raster layer.");  return;  }  // 展示图层  QgsProject::instance()->addMapLayer(rasterLayer);  }

5.4 栅格计算器

功能点描述：编写一个具有加、减、乘、除、对数、指数计算功能的栅格图层计算工具，计算结果以tiff格式的栅格文件输出。

操作及截图：

通过单击输入栅格图层按钮，调用on_InputpushButton_clicked()函数。函数中使用QFileDialog::getOpenFileName()函数打开文件对话框，让用户选择输入的TIFF文件。选择完成后，将文件路径显示在inputedit文本框中。然后使用GDAL库打开输入的TIFF图像文件，并将其作为GDAL数据集(GDALDataset)存储在inputDataset指针中。

之后通过组合框comboBox选择要计算的波段号。波段号从1开始，通过循环将波段号添加到组合框中。在计算表达式中输入进行计算的表达式；选择保存计算结果的位置：通过单击"OutputpushButton"按钮，调用on_OutputpushButton_clicked()函数。函数中使用QFileDialog::getSaveFileName()函数打开文件对话框，让用户选择保存计算结果的位置和文件名。选择完成后，将文件路径显示在outputedit文本框中。

最后点击"计算"按钮：通过单击"pushButton"按钮，调用on_pushButton_clicked()函数。函数中获取输入文件路径、输出文件路径和计算表达式。然后使用GDAL库打开输入的TIFF图像文件，并创建输出的TIFF图像文件。

图 96 栅格计算器演示界面

图 97 进行计算的源图层

图 98 栅格计算后的成果

核心代码：

void tiffcalculator::on_InputpushButton_clicked() {  QString filePath = QFileDialog::getOpenFileName(nullptr, "Select tiff", "", "Shapefiles (*.img *.tif *.tiff)");  ui.inputedit->setText(filePath);  GDALAllRegister();  // 注册GDAL驱动  // 打开输入数据集  GDALDataset* inputDataset = (GDALDataset*)GDALOpen(filePath.toUtf8().constData(), GA_ReadOnly);  if (inputDataset == nullptr) {  QMessageBox::critical(this, "错误", "无法打开输入文件。");  return;  }  // 使用波段号填充组合框  int bandCount = inputDataset->GetRasterCount();  ui.comboBox->clear();  for (int i = 1; i <= bandCount; ++i) {  ui.comboBox->addItem(QString::number(i));  }  GDALClose(inputDataset);  
}  void tiffcalculator::on_OutputpushButton_clicked() {  QString saveFilePath = QFileDialog::getSaveFileName(nullptr, "Save tiff file", "", "files (*.tiff)");  ui.outputedit->setText(saveFilePath);  
}  void tiffcalculator::on_pushButton_clicked() {  QString inputPath = ui.inputedit->text();  QString outputPath = ui.outputedit->text();  QString expression = ui.calculatelineEdit->text();  int selectedBand = ui.comboBox->currentText().toInt(); // 获取选择的波段号  GDALAllRegister();  // 注册GDAL驱动  // 打开输入图层  GDALDataset* inputDataset = (GDALDataset*)GDALOpen(inputPath.toUtf8().constData(), GA_ReadOnly);  if (inputDataset == nullptr) {  QMessageBox::critical(this, "错误", "打开输入图层失败");  return;  }  // 创建输出图层  GDALDriver* outputDriver = GetGDALDriverManager()->GetDriverByName("GTiff");  GDALDataset* outputDataset = outputDriver->Create(outputPath.toUtf8().constData(),inputDataset->GetRasterXSize(),inputDataset->GetRasterYSize(),1,GDT_Float32,nullptr);  if (outputDataset == nullptr) {  QMessageBox::critical(this, "错误", "创建输出图层失败");  GDALClose(inputDataset);  return;  }  // 获取输入图层的宽度和高度  int width = inputDataset->GetRasterXSize();  int height = inputDataset->GetRasterYSize();  // 读取输入图层数据  float* inputData = new float[width * height];  inputDataset->GetRasterBand(selectedBand)->RasterIO(GF_Read, 0, 0, width, height, inputData, width, height, GDT_Float32, 0, 0);  // 计算  float* outputData = new float[width * height];  for (int i = 0; i < width * height; ++i) {  float x = inputData[i];  float result = evaluateExpression(expression, x);  if (std::isnan(result)) {  QMessageBox::critical(this, "错误", "无效的表达式");  delete[] inputData;  delete[] outputData;  GDALClose(inputDataset);  GDALClose(outputDataset);  return;  }  outputData[i] = result;  }  // 将计算结果写入输出图层  outputDataset->GetRasterBand(1)->RasterIO(GF_Write, 0, 0, width, height, outputData, width, height, GDT_Float32, 0, 0);  // 释放内存并关闭文件  delete[] inputData;  delete[] outputData;  GDALClose(inputDataset);  GDALClose(outputDataset);  close();  QMessageBox::information(this, "成功", "计算完成！");  // 打开计算后的栅格文件  QFileInfo fi(outputPath);  if (!fi.exists()) { return; }  QString layerBaseName = fi.baseName(); // 图层名称  QgsRasterLayer* rasterLayer = new QgsRasterLayer(outputPath, layerBaseName);  if (!rasterLayer->isValid())  {  QMessageBox::critical(this, "错误", "打开文件失败");  return;  }  // 展示图层  QgsProject::instance()->addMapLayer(rasterLayer);  
}