基於GDAL庫，讀取海洋風場資料（.nc格式）c++版
阿新 • • 發佈：2018-12-03
經過這一段時間的對海洋資料的處理，接觸了大量的與海洋相關的資料，例如海洋地形、海洋表面溫度、鹽度、溼度、雲場、風場等資料，除了地形資料是grd格式外，其他的都是nc格式的資料。本文將以海洋風場資料為例，進行nc格式檔案的讀取。
海洋風場資料（ccmp_wind）一般情況下會包含三個資料集：第一個資料集是uwnd（standard_name = "eastward_wind"），第二個資料集是vwnd（standard_name = "northward_wind"），第三個資料集是nobs或者wspd。前兩個資料集是向量資料，表示此處的風場方向最後一個數據集是標量資料，代表此處的風速。每個資料集中資料的儲存又分為四個波段（也可以說是圖層），一天的觀測時間分為四個時間點，所以有四個圖層。
GDAL庫可以提供對nc格式資料的讀取，本次資料的讀取是在qt+vs2017環境下配置gdal庫和netcdf庫，環境的配置可以在網上找到，GDAL庫的配置可以根據《GDAL原始碼剖析和開發指南》書中的內容進行編譯和配置，配置完成後就可以執行資料，讀取nc檔案。
資料讀取的程式碼如下：
標頭檔案：
1 #ifndef CCMPFILEREAD_H 2 #define CCMPFILEREAD_H 3 class ccmpFileRead4 { 5 public: 6 void ccmpFileRead::fileread(const char*ccmpFilename); 7 }; 8 9 10 11 #endif // CCMPFILEREAD_H
原始檔：
1 #include "ccmpfileread.h" 2 3 #include <gdal_priv.h> 4 #include <vector> 5 #include <QVector> 6 7 #include <string> 8 #include <QString> 9#include <QStringList> 10 #include <QDebug> 11 12 #include <fstream> 13 14 using namespace std; 15 16 void ccmpFileRead::fileread(const char *ccmpFilename) 17 { 18 vector <string> vFileSets; 19 vector <string> pStrDesc; 20 vector<vector<float>> allSSTPixelNum1,allSSTPixelNum2,allSSTPixelNum3; 21 22 23 GDALAllRegister(); 24 CPLSetConfigOption("GDAL_FILENAME_IS_UTF8","NO");//中文路徑 25 GDALDataset* fileDataset = (GDALDataset*) GDALOpen(ccmpFilename,GA_ReadOnly);//開啟HDF資料集 26 if (fileDataset == NULL) 27 { 28 return; 29 } 30 31 char** sublist = GDALGetMetadata((GDALDatasetH) fileDataset,"SUBDATASETS");//獲得資料的字串，可以打印出來看看自己需要的資料在那 32 33 int iCount = CSLCount(sublist); 34 if(iCount <= 0){ 35 qDebug() << "該檔案沒有子資料" << endl; 36 GDALClose((GDALDriverH)fileDataset); 37 } 38 39 //儲存資料集資訊 40 for(int i = 0; sublist[i] != NULL;i++) 41 { 42 43 qDebug() << sublist[i] << endl; 44 45 if(i%2 != 0) 46 { 47 continue; 48 } 49 50 //三個資料集：uwnd vwnd wspd 只讀取前兩個資料集，第三個資料集是補充資料集 51 52 string tmpstr = sublist[i]; 53 tmpstr = tmpstr.substr(tmpstr.find_first_of("=")+1); 54 const char *tmpc_str = tmpstr.c_str(); 55 56 string tmpdsc = sublist[i+1]; 57 tmpdsc = tmpdsc.substr(tmpdsc.find_first_of("=")+1); 58 59 GDALDataset* hTmpDt = (GDALDataset*)GDALOpen(tmpc_str,GA_ReadOnly);//開啟該資料 60 61 if (hTmpDt != NULL) 62 { 63 vFileSets.push_back(tmpc_str); 64 } 65 if(&pStrDesc != NULL){ 66 pStrDesc.push_back(tmpdsc); 67 } 68 GDALClose(hTmpDt); 69 } 70 71 72 //三個資料集分別讀取 73 74 qDebug() << "read uwnd ......" << endl; 75 76 QString qtmpdsc1 = QString::fromStdString(pStrDesc[0]);//鎖定某一個數據集 77 78 qDebug()<<qtmpdsc1<<endl; 79 80 float *lineData = NULL; 81 if (qtmpdsc1!=NULL) 82 { 83 GDALDataset *tempDt = (GDALDataset *)GDALOpen(vFileSets[0].data(), GA_ReadOnly); 84 int BandNum = tempDt->GetRasterCount(); 85 86 int panBandmap[1] ={1}; 87 lineData = new float[1 * 200*200]; 88 tempDt->RasterIO(GF_Read,508,112,32,24,lineData,50,50,GDT_Float32,1,panBandmap,0,0,0); 89 90 91 for (int iLine = 0; iLine <tempDt->GetRasterYSize(); iLine++) 92 { 93 allSSTPixelNum1.resize(tempDt->GetRasterYSize()); 94 for (int iPixel = 0; iPixel < tempDt->GetRasterXSize(); iPixel++) 95 { 96 allSSTPixelNum1[iLine].resize(tempDt->GetRasterXSize()); 97 tempDt->RasterIO(GF_Read, 0, iLine, tempDt->GetRasterXSize(), 1,lineData, tempDt->GetRasterXSize(), 1, GDT_Float32, 1, panBandmap,0,0,0); 98 allSSTPixelNum1[iLine][iPixel] = lineData[iPixel]; 99 } 100 101 } 102 if(lineData) 103 { 104 delete[]lineData; 105 lineData = NULL; 106 } 107 108 qDebug() << "uwnd read over!" << endl; 109 110 qDebug() <<"uwnd="<<'\n'<<allSSTPixelNum1[200]<<'\n'<<endl; 111 112 } 113 114 //d讀取vwnd資料集 115 116 QString qtmpdsc2 = QString::fromStdString(pStrDesc[2]); 117 118 if (qtmpdsc2!=NULL) 119 { 120 GDALDataset *tempDt = (GDALDataset *)GDALOpen(vFileSets[0].data(), GA_ReadOnly); 121 int BandNum = tempDt->GetRasterCount(); 122 qDebug()<<BandNum<<endl; 123 int panBandmap[1] ={1}; 124 lineData = new float[1 * 200*200]; 125 tempDt->RasterIO(GF_Read,508,112,32,24,lineData,50,50,GDT_Float32,1,panBandmap,0,0,0); 126 127 128 for (int iLine = 0; iLine <tempDt->GetRasterYSize(); iLine++) 129 { 130 allSSTPixelNum2.resize(tempDt->GetRasterYSize()); 131 for (int iPixel = 0; iPixel < tempDt->GetRasterXSize(); iPixel++) 132 { 133 allSSTPixelNum2[iLine].resize(tempDt->GetRasterXSize()); 134 tempDt->RasterIO(GF_Read, 0, iLine, tempDt->GetRasterXSize(), 1,lineData, tempDt->GetRasterXSize(), 1, GDT_Float32, 1, panBandmap,0,0,0); 135 allSSTPixelNum2[iLine][iPixel] = lineData[iPixel]; 136 } 137 138 } 139 if(lineData) 140 { 141 delete[]lineData; 142 lineData = NULL; 143 } 144 145 qDebug() << "vwnd read over!" << endl; 146 147 qDebug() <<"vwnd="<<'\n'<<allSSTPixelNum2[200]<<'\n'<<endl; 148 149 } 150 151 //讀取wspd資料 152 153 QString qtmpdsc3 = QString::fromStdString(pStrDesc[2]); 154 155 if (qtmpdsc3!=NULL) 156 { 157 GDALDataset *tempDt = (GDALDataset *)GDALOpen(vFileSets[0].data(), GA_ReadOnly); 158 int BandNum = tempDt->GetRasterCount(); 159 qDebug()<<BandNum<<endl; 160 int panBandmap[1] ={1}; 161 lineData = new float[1 * 200*200]; 162 tempDt->RasterIO(GF_Read,508,112,32,24,lineData,50,50,GDT_Float32,1,panBandmap,0,0,0); 163 164 165 for (int iLine = 0; iLine <tempDt->GetRasterYSize(); iLine++) 166 { 167 allSSTPixelNum3.resize(tempDt->GetRasterYSize()); 168 for (int iPixel = 0; iPixel < tempDt->GetRasterXSize(); iPixel++) 169 { 170 allSSTPixelNum3[iLine].resize(tempDt->GetRasterXSize()); 171 tempDt->RasterIO(GF_Read, 0, iLine, tempDt->GetRasterXSize(), 1,lineData, tempDt->GetRasterXSize(), 1, GDT_Float32, 1, panBandmap,0,0,0); 172 allSSTPixelNum3[iLine][iPixel] = lineData[iPixel]; 173 } 174 175 } 176 177 if(lineData) 178 { 179 delete[]lineData; 180 lineData = NULL; 181 } 182 183 qDebug() << "wspd read over!" << endl; 184 185 qDebug() <<"wspd="<<'\n'<<allSSTPixelNum3[200]<<'\n'<<endl; 186 187 GDALClose((GDALDatasetH)tempDt); 188 189 } 190 191 GDALClose((GDALDriverH)fileDataset); 192 }
主函式呼叫：
1 #include <QCoreApplication> 2 #include <ccmpfileread.h> 3 int main(int argc, char *argv[]) 4 { 5 QCoreApplication a(argc, argv); 6 ccmpFileRead a1; 7 a1.fileread("E:/odp_workplace/odp_data/testdata/CCMP_Wind_Analysis_198707_V02.0_L3.5_RSS.nc"); 8 return a.exec(); 9 }
輸出結果：
如上圖所示資料已經讀取並顯示成功。