读取BIL格式高光谱数据——C/C++

在上一篇博客中，我们提到如何读取头文件。这篇博客将要谈到如何进一步读取高光谱数据本身，这里以BIL格式为例。
什么是BIL呢？BIL的全称为Band Interleave by Line，连续保存的是一个波段的一整行数据。如图1所示的三个波段的数据

图1

波段1的第一行为AAAAA，波段2的第一行为BBBBB，波段3的第一行为CCCCCC，那么整个第一行的数据就是AAAAABBBBBCCCCC。同理，第二行的数据为AAAAABBBBBCCCCC，第三行的数据为AAAAABBBBBCCCCC。再把3行的数据整体排开，得到AAAAABBBBBCCCCCAAAAABBBBBCCCCCAAAAABBBBBCCCCC，这也是写到二进制文件的顺序。
知道排列顺序是读取数据的关键，我们读取数据主要分为两个部分：一是读取空间位置点的光谱信息，二是读取某一波段的空间信息。为此，分别编写两个函数readSpec和readBand，分别返回一维指针和二维指针。读取的数据仍是上次的cupeff95.int，头文件的数据作为已知值给出。
具体代码如下所示：

// read_envi.cpp// C++ 读取ENVI的BIL 2-byte高光谱数据#include <iostream>using namespace std;int main(int argc, char const *argv[]){    printf("Hello, welcome to use C++!\n");    FILE *fp;    int lines = 350;     // 行数    int samples = 400;   // 列数    int bands = 50;      // 波段数    // 定义函数readSpec    short *readSpec(short *bil_data, int row, int col);     /* 返回的是指针 */    short *readBand(short *bil_data, int band_num);         /* 返回二维数组的指针 */    // 读取文件数据    fp = fopen("cup95eff.int", "rb");    int size = lines * samples * bands;    short *data = new short[size];    fread(data, sizeof(short), size, fp);    // 读取(i,j)的光谱数据    int i = 0, j = 0;    cout << "请输入读取的行、列：" << endl;    cin >> i >> j;    if (i > lines-1 || i < 0 || j > samples-1 || j < 0)    {        printf("超出范围\n");        exit(0);    }    short *sij = new short[50];    sij = readSpec(data, i, j);    // 输出光谱数据    printf("第%d行，第%d列的光谱数据：\n", i, j);    for (int i = 0; i < 50; ++i)    {        cout << sij[i] << endl;    }    // 读取指定波段    int k;    printf("请输入读取的波段号：\n");    cin >> k;    if (k > bands-1 || k < 0)    {        printf("超出范围\n");        exit(0);    }    short *band_k = new short[lines*samples];    band_k = readBand(data, k);    // 输出波段影像    printf("第%d波段影像：\n", k);    printf("输出左上角20x20区域：\n");    for (int i = 0; i < 20; ++i)    {        for (int j = 0; j < 20; ++j)        {            cout << *(band_k + i * 400 + j) << " ";        }        cout << endl;    }    return 0;}// readSpec函数：返回指定位置光谱数据short *readSpec(short *data, int i, int j){    short *spectrum = new short[50];    for (int k = 0; k < 50; ++k)    {        spectrum[k] = data[50 * 400 * i + j + 400 * k];    }    return spectrum;}// readBand函数：返回指定波段数据short *readBand(short *data, int k){    short *image = new short[350 * 400];    for (int i = 0; i < 350; ++i)    {        for (int j = 0; j < 400; ++j)        {            image[i * 400 + j] = data[400 * 50 * i + 400 * k + j];        }    }    return image;}

其中的关键就是把整个lines*samples*bands*2-byte大小的指针data进行分割，正如readSpec函数和readBand所做的一样。
编译、运行以上代码，读取并显示第99行，99列的光谱数据。读取第24波段的空间数据，显示左上角20x20区域的数据。
运行结果分别如图2、图3所示：

图2 读取光谱


图3 读取波段

我们用Matlab读取同样的数据来验证正确与否，matlab代码如下：

filename = 'cup95eff.int';rows = 350;cols = 400;bands = 50;data = multibandread(filename, [rows cols bands],...        'int16', 0, 'bil', 'ieee-le');Y = squeeze(data(:,:,25));Z = squeeze(data(100,200,:));

因为matlab的起始是1，C的起始是0，所以这里读的波段号为25，坐标为100。matlab的结果分别如图4、图5所示：

图4 光谱数据

图5 空间数据

对比可以得出，我们读取的数据正确。
读取BIL格式高光谱数据的工作基本完成！