百度360必应搜狗淘宝本站头条
当前位置:网站首页 > 技术分析 > 正文

C++借助GDAL批量读取栅格影像并生成像元数值的时间序列

liebian365 2025-01-26 23:08 14 浏览 0 评论

本文介绍基于C++ 语言GDAL库,批量读取大量栅格遥感影像文件,并生成各像元数值的时间序列数组的方法。

首先,我们来明确一下本文所需实现的需求。现在有一个文件夹,其中包含了很多不同格式的文件,如下图所示。

其中,我们首先需要遍历这一文件夹,遴选出其中所有类型为.bmp格式的栅格遥感影像文件(一共有6个),并分别读取文件(已知这些遥感影像的行数、列数都是一致的);随后,将不同遥感影像同一个位置的像素的数值进行分别读取,并存储在一个数组中。例如,最终我们生成的第一个数组,其中共有6个元素,分别就是上图所示文件夹中6景遥感影像各自(0,0)位置的像元数值;生成的第二个数组,其中也是6个元素,分别就是6景遥感影像各自(1,0)位置的像元数值,以此类推。其中,显然我们得到的数组个数,就是遥感影像像元的个数。此外,这里6景遥感影像的排序,是按照文件名称的升序来进行的。

明确了具体需求,接下来就可以开始代码的实践。其中,本文分为两部分,第一部分为代码的分段讲解,第二部分为完整代码。

此外,本文是基于GDAL库来实现栅格数据读取的;具体GDAL库的配置方法大家可以参考文章Visual Studio配置C++中GDAL、SQLite与PROJ环境并编译的方法

1 代码分段介绍

1.1 代码准备

这一部分主要是代码的头文件命名空间与我们自行撰写的自定义函数get_need_file()的声明;具体代码如下所示。

#include <iostream>
#include <vector>
#include <io.h>
#include "gdal_priv.h"

using namespace std;

void get_need_file(string path, vector<string>& file, string ext);

其中,由于我们在接下来的代码中需要用到容器vector这一数据类型,因此首先需要添加#include <vector>;同时,我们在接下来的代码中需要用到头文件io.h中的部分函数(主要都是一些与计算机系统、文件管理相关的函数),因此需要添加#include <io.h>;此外,我们是基于GDAL库来实现栅格数据读取的,因此需要添加#include "gdal_priv.h"

接下来,这里声明了一个自定义函数get_need_file(),具体我们在本文1.2部分介绍。

1.2 栅格文件筛选

由于我这里几乎将全部的代码都放在了主函数中,因此这一部分就先介绍代码main()函数的第一部分,亦即栅格文件的遴选部分;具体代码如下所示。

int main() {
    string file_path = R"(E:\02_Project\02_ChlorophyllProduce\01_Data\00_Test)";
    vector<string> my_file;
    string need_extension = ".bmp";
    get_need_file(file_path, my_file, need_extension);
    int file_size = my_file.size();
    if (file_size == 0)
    {
        cout << "No file can be found!" << endl;
    }
    else
    {
        cout << "Find " << file_size << " file(s).\n" << endl;
    }

这一部分主要就是做好调用自定义函数get_need_file()的变量准备,并调用get_need_file()函数,得到指定文件夹下的栅格文件;随后,将栅格文件的筛选结果进行输出。这一部分的具体代码介绍,大家查看文章C++遍历文件夹筛选出指定格式的文件或具有特定名称的文件即可,这里就不再赘述。

1.3 栅格文件读取

这一部分主要是基于GDAL库,循环读取前述文件夹中的每一个栅格遥感影像文件。

    int nXSize, nYSize;
    float** pafScanline = new float* [file_size];
    int pic_index = 1;
    for (auto x : my_file)
    {
        GDALDataset* poDataset;
        GDALAllRegister();
        CPLSetConfigOption("GDAL_FILENAME_IS_UTF8", "NO");
        poDataset = (GDALDataset*)GDALOpen(x.c_str(), GA_ReadOnly);
        if (poDataset == NULL)
        {
            cout << "Open File " << x << " Error!" << endl;
        }
        else
        {
            cout << "Open File " << x << " Success!" << endl;
        }

        GDALRasterBand* poBand;
        poBand = poDataset->GetRasterBand(1);
        nXSize = poBand->GetXSize();
        nYSize = poBand->GetYSize();
        cout << nXSize << "," << nYSize << "\n" << endl;
        pafScanline[pic_index - 1] = new float[nXSize * nYSize];
        poBand->RasterIO(GF_Read, 0, 0, nXSize, nYSize, pafScanline[pic_index - 1], nXSize, nYSize, GDT_Float32, 0, 0);

        pic_index ++;
    }

其中,nXSizenYSize分别表示栅格遥感影像的列数与行数,pafScanline是我们读取栅格遥感影像文件所需的变量,之后读取好的遥感影像数据就会存放在这里;由于我们有多个栅格文件需要读取,因此通过for循环来实现批量读取的操作,并通过pic_index这个变量作为每一次读取文件的计数。

在这里,float** pafScanline = new float* [file_size];这句代码表示我们将pafScanline作为一个指向指针的指针的数组;在后期读取遥感影像数据后,pafScanline[0]pafScanline[1]一直到pafScanline[5],这6个数值同样分别是指针,分别指向存储6景遥感影像数据的地址。这里我们通过new实现对pafScanline内存的动态分配,因为我们在获取栅格遥感影像的景数(也就是文件夹中栅格遥感影像文件的个数)之前,也不知道具体需要给pafScanline这一变量分配多少的内存。此外,在for循环中,我们还对pafScanline[0]pafScanline[1]一直到pafScanline[5]同样进行了动态内存分配,因为我们在获取每一景栅格遥感影像的行数与列数之前,同样是不知道需要给pafScanline[x]6个数组变量分配多少内存的。

随后,for循环中的其他部分,就是GDAL库读取遥感影像的基本代码。读取第一景遥感影像数据后,我们将数据保存至pafScanline[0],并随后进行第二次循环,读取第二景遥感影像数据,并将其数据保存至pafScanline[1]中,随后再次循环;以此类推,直至读取6景遥感影像完毕。

如果大家只是需要实现C++ 批量读取栅格遥感影像数据,那么以上操作就已经实现了大家的需求。其中,显然pafScanline[0]就是第一景遥感影像数据,pafScanline[1]就是第二景遥感影像数据,pafScanline[2]就是第三景遥感影像数据,以此类推。

1.4 像元时间序列数组生成

这一部分则是基于以上获取的各景遥感影像数据读取结果,进行每一个像元数值的时间序列数组生成。

    float** pixel_paf = new float* [nXSize * nYSize];
    for (int pixel_num = 0; pixel_num < nXSize * nYSize; pixel_num++)
    {
        pixel_paf[pixel_num] = new float[file_size];
        for (int time_num = 0; time_num < file_size; time_num++)
        {
            pixel_paf[pixel_num][time_num] = pafScanline[time_num][pixel_num];
        }
    }

这一部分的代码思路其实也非常简单,就是通过两个for循环,将原本一共6的、每一个表示每一景遥感影像中全部数据的数组,转变为一共X的(X表示每一景遥感影像的像元总个数)、每一个表示每一个位置的像元在6景遥感影像中的各自数值的数组。

在这里,由于同样的原因,我们对pixel_paf亦进行了内存的动态分配。

1.5 输出测试与代码收尾

这一部分主要是输出一个我们刚刚配置好的像元数值时间序列数组,从而检查代码运行结果是否符合我们的要求;此外,由于前面我们对很多变量进行了动态内存分配,因此需要将其delete掉;同时,这里还可以对前面我们定义的指向指针的指针赋值为NULL,这样子其就不能再指向任何地址了,即彻底将其废除。

    for (int i = 0; i < file_size; i++)
    {
        cout << pixel_paf[0][i] << "," << endl;
    }

    delete[] pafScanline;
    delete[] pixel_paf;
    pafScanline = NULL;
    pixel_paf = NULL;

    return 0;
}

至此,代码的主函数部分结束。

1.6 自定义函数

这一部分是我们的自定义函数get_need_file()

void get_need_file(string path, vector<string>& file, string ext)
{
    intptr_t file_handle = 0;
    struct _finddata_t file_info;
    string temp;
    if ((file_handle = _findfirst(temp.assign(path).append("/*" + ext).c_str(), &file_info)) != -1)
    {
        do
        {
            file.push_back(temp.assign(path).append("/").append(file_info.name));
        } while (_findnext(file_handle, &file_info) == 0);
        _findclose(file_handle);
    }
}

如前所述,这一部分的具体代码介绍,大家查看文章C++遍历文件夹筛选出指定格式的文件或具有特定名称的文件即可,这里就不再赘述。

2 完整代码

本文所需用到的完整代码如下所示。

#include <iostream>
#include <vector>
#include <io.h>
#include "gdal_priv.h"

using namespace std;

void get_need_file(string path, vector<string>& file, string ext);

int main() {
    string file_path = R"(E:\02_Project\02_ChlorophyllProduce\01_Data\00_Test)";
    vector<string> my_file;
    string need_extension = ".bmp";
    get_need_file(file_path, my_file, need_extension);
    int file_size = my_file.size();
    if (file_size == 0)
    {
        cout << "No file can be found!" << endl;
    }
    else
    {
        cout << "Find " << file_size << " file(s).\n" << endl;
    }

    int nXSize, nYSize;
    float** pafScanline = new float* [file_size];
    int pic_index = 1;
    for (auto x : my_file)
    {
        GDALDataset* poDataset;
        GDALAllRegister();
        CPLSetConfigOption("GDAL_FILENAME_IS_UTF8", "NO");
        poDataset = (GDALDataset*)GDALOpen(x.c_str(), GA_ReadOnly);
        if (poDataset == NULL)
        {
            cout << "Open File " << x << " Error!" << endl;
        }
        else
        {
            cout << "Open File " << x << " Success!" << endl;
        }

        GDALRasterBand* poBand;
        poBand = poDataset->GetRasterBand(1);
        nXSize = poBand->GetXSize();
        nYSize = poBand->GetYSize();
        cout << nXSize << "," << nYSize << "\n" << endl;
        pafScanline[pic_index - 1] = new float[nXSize * nYSize];
        poBand->RasterIO(GF_Read, 0, 0, nXSize, nYSize, pafScanline[pic_index - 1], nXSize, nYSize, GDT_Float32, 0, 0);

        pic_index ++;
    }

    float** pixel_paf = new float* [nXSize * nYSize];
    for (int pixel_num = 0; pixel_num < nXSize * nYSize; pixel_num++)
    {
        pixel_paf[pixel_num] = new float[file_size];
        for (int time_num = 0; time_num < file_size; time_num++)
        {
            pixel_paf[pixel_num][time_num] = pafScanline[time_num][pixel_num];
        }
    }

    for (int i = 0; i < file_size; i++)
    {
        cout << pixel_paf[0][i] << "," << endl;
    }

    delete[] pafScanline;
    delete[] pixel_paf;
    pafScanline = NULL;
    pixel_paf = NULL;

    return 0;
}

void get_need_file(string path, vector<string>& file, string ext)
{
    intptr_t file_handle = 0;
    struct _finddata_t file_info;
    string temp;
    if ((file_handle = _findfirst(temp.assign(path).append("/*" + ext).c_str(), &file_info)) != -1)
    {
        do
        {
            file.push_back(temp.assign(path).append("/").append(file_info.name));
        } while (_findnext(file_handle, &file_info) == 0);
        _findclose(file_handle);
    }
}

当我们运行上述代码后,将会出现如下所示的界面。

其中,会显示栅格遥感影像文件的筛选情况、具体文件名称及其各自的行号与列号;同时,最后一部分则是本文1.5部分提及的测试输出结果,其表示本文所用的6景遥感影像各自(0,0)位置处的像元数值。

至此,大功告成。

相关推荐

4万多吨豪华游轮遇险 竟是因为这个原因……

(观察者网讯)4.7万吨豪华游轮搁浅,竟是因为油量太低?据观察者网此前报道,挪威游轮“维京天空”号上周六(23日)在挪威近海发生引擎故障搁浅。船上载有1300多人,其中28人受伤住院。经过数天的调...

“菜鸟黑客”必用兵器之“渗透测试篇二”

"菜鸟黑客"必用兵器之"渗透测试篇二"上篇文章主要针对伙伴们对"渗透测试"应该如何学习?"渗透测试"的基本流程?本篇文章继续上次的分享,接着介绍一下黑客们常用的渗透测试工具有哪些?以及用实验环境让大家...

科幻春晚丨《震动羽翼说“Hello”》两万年星间飞行,探测器对地球的最终告白

作者|藤井太洋译者|祝力新【编者按】2021年科幻春晚的最后一篇小说,来自大家喜爱的日本科幻作家藤井太洋。小说将视角放在一颗太空探测器上,延续了他一贯的浪漫风格。...

麦子陪你做作业(二):KEGG通路数据库的正确打开姿势

作者:麦子KEGG是通路数据库中最庞大的,涵盖基因组网络信息,主要注释基因的功能和调控关系。当我们选到了合适的候选分子,单变量研究也已做完,接着研究机制的时便可使用到它。你需要了解你的分子目前已有哪些...

知存科技王绍迪:突破存储墙瓶颈,详解存算一体架构优势

智东西(公众号:zhidxcom)编辑|韦世玮智东西6月5日消息,近日,在落幕不久的GTIC2021嵌入式AI创新峰会上,知存科技CEO王绍迪博士以《存算一体AI芯片:AIoT设备的算力新选择》...

每日新闻播报(September 14)_每日新闻播报英文

AnOscarstatuestandscoveredwithplasticduringpreparationsleadinguptothe87thAcademyAward...

香港新巴城巴开放实时到站数据 供科技界研发使用

中新网3月22日电据香港《明报》报道,香港特区政府致力推动智慧城市,鼓励公私营机构开放数据,以便科技界研发使用。香港运输署21日与新巴及城巴(两巴)公司签署谅解备忘录,两巴将于2019年第3季度,开...

5款不容错过的APP: Red Bull Alert,Flipagram,WifiMapper

本周有不少非常出色的app推出,鸵鸟电台做了一个小合集。亮相本周榜单的有WifiMapper's安卓版的app,其中包含了RedBull的一款新型闹钟,还有一款可爱的怪物主题益智游戏。一起来看看我...

Qt动画效果展示_qt显示图片

今天在这篇博文中,主要实践Qt动画,做一个实例来讲解Qt动画使用,其界面如下图所示(由于没有录制为gif动画图片,所以请各位下载查看效果):该程序使用应用程序单窗口,主窗口继承于QMainWindow...

如何从0到1设计实现一门自己的脚本语言

作者:dong...

三年级语文上册 仿写句子 需要的直接下载打印吧

描写秋天的好句好段1.秋天来了,山野变成了美丽的图画。苹果露出红红的脸庞,梨树挂起金黄的灯笼,高粱举起了燃烧的火把。大雁在天空一会儿写“人”字,一会儿写“一”字。2.花园里,菊花争奇斗艳,红的似火,粉...

C++|那些一看就很简洁、优雅、经典的小代码段

目录0等概率随机洗牌:1大小写转换2字符串复制...

二年级上册语文必考句子仿写,家长打印,孩子照着练

二年级上册语文必考句子仿写,家长打印,孩子照着练。具体如下:...

一年级语文上 句子专项练习(可打印)

...

亲自上阵!C++ 大佬深度“剧透”:C++26 将如何在代码生成上对抗 Rust?

...

取消回复欢迎 发表评论: