C++ GDAL找出多时相遥感影像缺失的日期并自动生成新的全零图像作为替补

简介: C++ GDAL找出多时相遥感影像缺失的日期并自动生成新的全零图像作为替补

  本文介绍基于C++ 语言的GDAL库,基于一个存储大量遥感影像文件夹,依据每一景遥感影像的文件名中表示日期的那个字段,找出这些遥感影像中缺失的成像日期,并新生成多个像元值全部为0的栅格文件,作为这些缺失日期当日的遥感影像文件的方法。

  首先,我们来看一下本文需要实现的需求。现在有一个文件夹,存储了从2018年第001天到2022年第361天的全部遥感影像,其中每一景图像的像元个数、空间参考信息、NoData值等都是一致的。对于这些遥感影像,原本应该是每10天就有1景;但是由于遥感影像数据有缺失,因此部分日期没有对应的遥感影像。如下图所示,可以看到比如2018年的061这一天,它就没有对应的遥感影像。

  但是,由于后期处理的需要,我们现在希望对这些缺失日期的遥感影像文件加以填补——具体的需求是,我们新建若干个像元值全部为0的栅格文件,作为每一个缺失日期当日的遥感影像文件;这些填补的、新的遥感影像文件的各项信息(比如像元个数、空间参考信息等)都和原本的文件一致即可,只要保证全部的像元都是0就行。

  知道了需求,我们就可以开始代码的撰写。本文用到的代码具体如下所示。其中,关于C++ 语言配置GDAL库的方法,大家可以参考文章Visual Studio配置并编译C++环境下GDAL库、SQLite环境与PROJ库的方法https://blog.csdn.net/zhebushibiaoshifu/article/details/127088090)。

#include <iostream>
#include <iomanip>
#include <sstream>
#include "gdal_priv.h"
#include "cpl_conv.h"
using namespace std;
void create_missing_raster(string path);
int main() {
    string file_path = R"(E:\02_Project\TIFF\TEST)";
    create_missing_raster(file_path);
    return 0;
}
void create_missing_raster(string path)
{
  vector<string> all_file_path;
  for (int year = 2018; year <= 2022; year++)
  {
    for (int day = 1; day <= 361; day += 10)
    {
      ostringstream oss;
      oss << path << "/Albers_MuSyQ.NDVI.16m." << to_string(year) << setfill('0') << setw(3) << to_string(day) << "000000.49SDB.001.h5NDVI-NDVI.tif";
      all_file_path.push_back(oss.str());
    }
  }
  GDALAllRegister();
  int x_size, y_size;
  string one_actual_path;
  for (string& one_file_path : all_file_path)
  {
    GDALDataset* poDataset_actual;
    poDataset_actual = (GDALDataset*)GDALOpen(one_file_path.c_str(), GA_ReadOnly);
    if (poDataset_actual != NULL)
    {
      one_actual_path = one_file_path;
      x_size = poDataset_actual->GetRasterXSize();
      y_size = poDataset_actual->GetRasterYSize();
      GDALClose((GDALDatasetH)poDataset_actual);
      break;
    }
  }
  for (string& one_file_path : all_file_path)
  {
    if (CPLCheckForFile((char*)one_file_path.c_str(), NULL) == FALSE)
    {
      GDALDataset* poDataset;
      GDALDriver* poDriver;
      poDriver = GetGDALDriverManager()->GetDriverByName("GTiff");
      GDALDataset* poDataset_actual;
      poDataset_actual = (GDALDataset*)GDALOpen(one_actual_path.c_str(), GA_ReadOnly);
      poDataset = poDriver->CreateCopy(one_file_path.c_str(), poDataset_actual, FALSE, NULL, NULL, NULL);
      double* pafScanline = new double[x_size * y_size];
      for (int i = 0; i < x_size * y_size; i++)
      {
        pafScanline[i] = 0.0;
      }
      GDALRasterBand* poBand;
      poBand = poDataset->GetRasterBand(1);
      poBand->RasterIO(GF_Write, 0, 0, x_size, y_size, pafScanline, x_size, y_size, GDT_Float64, 0, 0);
      delete[] pafScanline;
      GDALClose((GDALDatasetH)poDataset);
      GDALClose((GDALDatasetH)poDataset_actual);
      cout << "New file is :" << one_file_path << endl;
    }
  }
  GDALDestroyDriverManager();
}

  上述代码主要都是在create_missing_raster(string path)函数内实现具体功能的,我们主要就对这一函数加以讲解。

  首先,我们需要基于文件夹中遥感影像文件的文件名称特征,遍历生成文件名列表。在这里,我们使用两个嵌套的for循环,生成所有可能的栅格图像文件名,并将这些文件名保存在all_file_path向量中。其中,栅格图像的文件名根据年份和天数生成,并通过setfill('0')setw(3)这两个函数保证我们生成的日期满足YYYYDDD这种格式。

  随后,基于GDALAllRegister这一GDAL库的初始化函数,用于注册所有支持的数据格式驱动程序。接下来,我们使用GDALOpen函数,从2018001这一天开始,通过循环打开对应名字的文件,直到找到文件夹中第一个实际存在的栅格图像文件poDataset_actual),并获取其栅格图像的行列数(x_sizey_size);我们后期的操作需要用到这个行列数,并且会将这个实际存在的栅格文件作为生成新的栅格文件的模板。

  接下来,我们遍历文件名列表all_file_path,对每个文件名进行处理。对于不存在的栅格图像文件,使用GDALDriver创建一个新的数据集(poDataset),并将其中的像元值设置为0。如果栅格图像文件已经存在,则跳过不处理。其中,在对缺失的栅格图像加以生成时,我们首先使用GetGDALDriverManager()->GetDriverByName函数获取GDAL驱动程序对象,然后使用CreateCopy函数创建新的栅格图像;其中,我们就是以前期找到的文件夹中第一个实际存在的栅格图像文one_actual_path为模板。随后,我们用0填充新创建的栅格图像,并使用RasterIO函数对栅格图像的像元进行写入操作。

  最后,在上述处理完成后,使用GDALClose函数关闭数据集,并输出新创建的栅格图像的文件名。随后,我们使用GDALDestroyDriverManager销毁GDAL驱动程序管理器,释放资源。

  随后,我们运行代码,可以看到每一个新生成的栅格图像文件(也就是原本当日没有成像的日期对应的遥感影像)都会打印出来。

  随后,我们打开文件夹,可以看到之前没有遥感影像的日期,目前也都存在一景遥感影像与其对应了。比如2018年的061这一天,目前已经有了一景遥感影像。

  至此,大功告成。

欢迎关注:疯狂学习GIS

相关文章
|
1月前
|
编译器 C++
C++初阶--类与对象--const成员和日期类的实现
C++初阶--类与对象--const成员和日期类的实现
|
1月前
|
存储 C++
【C++初阶】第三站:类和对象(中) -- 日期计算器
【C++初阶】第三站:类和对象(中) -- 日期计算器
|
1月前
|
存储 C语言 C++
c++日期和时间
c++日期和时间
19 0
|
2月前
|
存储 数据采集 数据可视化
【C++】医院PACS医学图像存储和传输系统源码
图像后处理与重建 •MPR\CPR(三维多平面重建) •VRT(三维容积重建) •SSD(三维表面重建) •VE(虚拟内窥镜) •MIP(最大密度投影)、MinIP(最小密度投影) •CalSCore(心脏图像冠脉钙化积分)
31 3
|
28天前
|
算法 C语言 C++
C++ std::chrono库使用指南 (实现C++ 获取日期,时间戳,计时等功能)(一)
C++ std::chrono库使用指南 (实现C++ 获取日期,时间戳,计时等功能)
112 1
|
3月前
|
存储 数据处理 开发工具
Baumer工业相机堡盟工业相机如何通过NEOAPI SDK实现相机的高速图像保存(C++)
Baumer工业相机堡盟工业相机如何通过NEOAPI SDK实现相机的高速图像保存(C++)
41 0
|
28天前
|
缓存 算法 安全
C++ std::chrono库使用指南 (实现C++ 获取日期,时间戳,计时等功能)(三)
C++ std::chrono库使用指南 (实现C++ 获取日期,时间戳,计时等功能)
36 1
|
3月前
|
存储 数据处理 开发工具
Baumer工业相机堡盟工业相机如何通过NEOAPI SDK实现相机图像转换由Mono10转换为Mono8(C++)
Baumer工业相机堡盟工业相机如何通过NEOAPI SDK实现相机图像转换由Mono10转换为Mono8(C++)
45 0
|
28天前
|
安全 调度 C++
C++ std::chrono库使用指南 (实现C++ 获取日期,时间戳,计时等功能)(二)
C++ std::chrono库使用指南 (实现C++ 获取日期,时间戳,计时等功能)
54 0
|
3月前
|
监控 API 开发工具
Baumer工业相机堡盟工业相机如何通过NEOAPI SDK获取每张图像的微秒时间和FrameID功能(C++)
Baumer工业相机堡盟工业相机如何通过NEOAPI SDK获取每张图像的微秒时间和FrameID功能(C++)
33 0