遥感处理与应用

C++读取NC格式文件的注意事项

2024-06-18  本文已影响0人  疯狂学习GIS

  本文介绍基于C++语言的netCDF库读取.nc格式的栅格文件时,出现数据无法读取数据读取错误无法依据维度提取变量等情况的原因与解决方法。

  最近,由于需要读取ERA5气象数据,因此使用C++语言中的netCDF库读取.nc格式文件;这其中也是踩了很多的坑,所以在此记录一下,也相当于是汇总了C++netCDF库读取.nc文件时的一些常见问题。

1 环境配置

  环境配置是使用每一个第三方库时,首先遇到的问题。针对不同语言,netCDF库也有着对应的不同版本,我们之前的文章也介绍过在PythonC++等代码的开发环境中,配置netCDF库的具体方法,大家直接参考文章Python模块h5py与netCDF4在Anaconda环境的下载与安装,以及文章部署C++中netCDF库:Visual Studio即可。但当然,本文后续提到的各类问题,都是在基于C++语言的netCDF库读取.nc格式文件时,出现的问题(毕竟Python读取.nc格式文件还是很方便的,感觉一般也不会有太大的问题或坑)。

2 数据增益值与偏移值

  首先,如果大家读取数据时,发现得到的结果数值很奇怪、不符合数据范围的实际情况,那么大概率就是.nc文件的变量存在scaleoffset(增益值、偏移值)导致的;关于这一点,我们之前已经用了完整的一篇文章对其加以介绍,大家参考文章C++提取NC文件时结果错误的一种解决思路即可。

3 NcVar格式数据不能跨函数使用

  此外,经过操作发现,我们读取得到的NcVar格式的变量,其似乎只能在.nc格式文件被读取后立刻使用,而不能跨函数使用;如果跨函数使用,就会出现badid的错误。

  例如,我们一般情况下,都是通过如下代码,打开.nc格式文件,并读取其中的变量数据。

NcFile file(path, NcFile::read);
NcVar var = file.getVar("ssrd");

  那么此时,假设我上述的2句代码是在一个函数中运行的,然后我将得到的NcVar格式的变量var作为这个函数返回值,返回给调用者;随后,又将这个被返回的NcVar格式变量var作为参数,输入到另一个新的函数中——那么,在这个新的函数中,我们如果还想读取var这个变量(例如用如下的代码来读取var这个变量),就会出现badid的错误。

vector<size_t> start{ size_t(time_idx), size_t(latitude_idx), size_t(longitude_idx) };
double* value = new double[1];
var.getVar(start, value);

  在这个地方,一开始我还以为是我的start参数设置有误,导致一直无法读取var;后来才注意到,原来是这个var不能够跨越函数读取导致的(应该是这样的吧,具体我倒也没有看到官网上有明确的说明)。

4 时间维度需要放在第一个位置

  有时,我们需要按照不同维度,对变量数据加以读取。例如,我这里的.nc格式文件中,变量是1种气象数据,其具有3维度,包括经度、纬度与时间等。

  那么,假设我们希望获取某一个指定时间中,某一个经度与纬度处,对应的变量的数值(相当于就是其在栅格文件中的像素值)——那么多数情况下,我们会选择.getVar(start, value)这种方法,对变量数据加以读取。如果是如此,就需要注意将时间维度放在start第一个元素的位置上;具体代码如下所示。

NcFile file(path, NcFile::read);
NcVar var = file.getVar("ssrd");
vector<size_t> start{ size_t(time_idx), size_t(latitude_idx), size_t(longitude_idx) };
double* value = new double[1];
var.getVar(start, value);

  可以看到,我们需要保证上述代码中,time_idx是这个start变量的第一个元素。关于这一点,在C++版本的netCDF库的官方网站中,也有具体提及,如下图所示。

  在这里多提一句。我们可以用下述代码,将读取.nc格式的栅格文件,并获取其中的所有维度。这个代码首先读取了我们的.nc格式文件,然后通过.getDims()方法,获取了其中的全部维度,并将每一个维度都放入了multimap格式的变量dimension_map中。

NcFile file(path, NcFile::read);
multimap<string, NcDim> dimension_map;
dimension_map = file.getDims();

  此时,当我们查看这个dimension_map,有时会发现——时间维度并不在这个multimap的第一位(multimap是有序的),其myId参数的值也不一定是0;如下图所示。

  如下图所示,我打开了另一个.nc格式的栅格文件并查看其维度,可以看到虽然此时时间维度myId值为0,但是其在multimap中的位置依然不是第一位。

  但是,尽管如此——只要我们需要基于.getVar(start, value)这种方法,对变量数据加以读取,那么就一定注意将时间维度放在start的第一个位置。

  至此,大功告成。

上一篇下一篇

猜你喜欢

热点阅读