1 前言
纹理特征在地物光谱特征比较相似的时候常作为一种特征用于图像的分类和信息提取。本文主要介绍基于灰度共生矩阵的图像纹理提取。
2 灰度共生矩阵
纹理石油灰度分布在空间位置上反复出现而形成的，因而图像空间中相隔某距离的两个像素之间存在一定的灰度关系，即图像中灰度的空间相关特性。灰度共生矩阵是一种通过研究灰度的空间相关特性来描述纹理的常用方法。
有一个网站提供了非常详细的关于灰度共生矩阵的介绍，我在这里就不赘述了。
3 纹理提取的算法流程
本文重点说明一下遥感影像纹理提取的算法流程，如下图1所示：
图 1 纹理提取算法流程
具体描述如下：
1) 灰度降级，对原始影像进行灰度降级如8,16,32,64等；
纹理计算的灰度降级策略来源于IDL的bytscl函数介绍，具体描述如下：
图 2 灰度降级
2) 根据设定好的窗口大小，逐窗口计算灰度共生矩阵；
3) 根据选择的二阶统计量，计算纹理值。
4 纹理算子
协同性（GLCM_HOM）：对应ENVI的Homogeneity
反差性（GLCM_CON）：
非相似性（GLCM_DIS）：
均值GLCM_MEAN：对应ENVI的Mean
方差GLCM_VAR：对应ENVI的Variance
角二阶矩GLCM_ASM：对应ENVI的Second Moment
相关性GLCM_COR：对应ENVI的Correlation
GLDV角二阶矩GLDV_ASM：
熵GLCM_ENTROPY：对应ENVI的Entropy
归一化灰度矢量均值GLDV_MEAN：对应ENVI的Dissimilarity
归一化角二阶矩GLDV_CON：对应ENVI的Contrast
5 代码实现
要是有哪位大神可以指点下如何提高纹理计算的效率感激不尽~
先上传头文件和源文件
1 #pragma once
2 #include <iostream>
3 #include <string>
4 #include <vector>
5
6 class CCo_occurrenceTextureExtration
7 {
8 public:
9 CCo_occurrenceTextureExtration(void);
10 ~CCo_occurrenceTextureExtration(void);
11
12
13 typedef enum ANGLE_TYPE
14 {
15 DEGREE_0,
16 DEGREE_45,
17 DEGREE_90,
18 DEGREE_135
19 };
20 typedef enum LIST_TYPE
21 {
22 GLCM_HOM,
23 GLCM_CON,
24 GLCM_DIS,
25 GLCM_MEAN,
26 GLCM_VAR,
27 GLCM_ASM,
28 GLCM_COR,
29 GLDV_ASM,
30 GLCM_ENTROPY,
31 GLDV_MEAN,
32 GLDV_CON,
33 };
34
35
36 //影像灰度降级（包括直方图均衡化）
37 bool Init();
38 bool ReduceGrayLevel();
39 bool ExtraTextureWinEXE();
40 bool ExtraTexturePixEXE();
41 int ExtraTexturePixOpenCLEXE();
42 private:
43 char* LoadProgSource(const char* cFilename, const char* cPreamble, size_t* szFinalLength);
44 bool GetGrayCoocurrenceMatrix(std::vector<std::vector<int>>grayValue,
45 std::vector<std::vector<float>> &grayMatrix,int &m_count);
46 float CalGLCMStatistics(std::vector<std::vector<float>> grayMatrix,
47 const int count);
48 double* cdf(int *h,int length);
49
50 public:
51 std::string m_strInFileName;
52 std::string m_strOutFileName;
53 std::string m_strReduceLevelFileName;
54 int m_iWidth;
55 int m_iHeigth;
56 int m_iBandCount;
57 int *m_BandMap;
58 int m_AngleMode;
59 int m_TextureMode;
60 int m_dis;
61 int m_grayLevel;
62 int m_winWidth;
63 int m_winHeigth;
64
65 };源文件如下：
1 #include "Co_occurrenceTextureExtration.h"
2 #include <gdal_alg_priv.h>
3 #include <gdal_priv.h>
4 #include <math.h>
5 #include <omp.h>
6 #include <CL/cl.h>
7 #include <stdlib.h>
8 #include <malloc.h>
9 #pragma comment(lib,"OpenCL.lib")
10
11 CCo_occurrenceTextureExtration::CCo_occurrenceTextureExtration( )
12 {
13 m_BandMap = NULL;
14 }
15 CCo_occurrenceTextureExtration::~CCo_occurrenceTextureExtration(void)
16 {
17 if(!m_BandMap){
18 delete []m_BandMap;
19 }
20 }
21
22 bool CCo_occurrenceTextureExtration::Init()
23 {
24 GDALAllRegister();
25 GDALDataset *poDS = (GDALDataset*)GDALOpen(m_strInFileName.c_str(),GA_ReadOnly);
26 m_iWidth = poDS->GetRasterXSize();
27 m_iHeigth = poDS->GetRasterYSize();
28 m_iBandCount = poDS->GetRasterCount();
29
30 m_BandMap = new int[m_iBandCount];
31 int i = 0;
32 while(i < m_iBandCount){
33 m_BandMap[i] = i+1;
34 i++;
35 }
36 GDALClose((GDALDatasetH)poDS);
37 return TRUE;
38 }
39
40 bool CCo_occurrenceTextureExtration::ReduceGrayLevel()
41 {
42 //直方图均衡化
43 ///// 统计直方图
44 GDALAllRegister();
45 GDALDataset *poDS = (GDALDataset*)GDALOpen(m_strInFileName.c_str(),GA_ReadOnly);
46 double adfMSSGeoTransform[6] = {0};
47 poDS->GetGeoTransform(adfMSSGeoTransform);
48 GDALDataType eDT = poDS->GetRasterBand(1)->GetRasterDataType();
49 //创建文件
50 GDALDriver *poOutDriver = (GDALDriver*)GDALGetDriverByName("GTIFF");
51 GDALDataset *poOutDS = poOutDriver->Create(m_strReduceLevelFileName.c_str(),
52 m_iWidth,m_iHeigth,m_iBandCount,GDT_Byte,NULL);
53 poOutDS->SetGeoTransform(adfMSSGeoTransform);
54 poOutDS->SetProjection(poDS->GetProjectionRef());
55 for(int i = 0;i<m_iBandCount;i++){
56 double *bandMINMAX = new double[2];
57 GDALRasterBand *poBand = poDS->GetRasterBand(m_BandMap[i]);
58 GDALRasterBand *poNewBand = poOutDS->GetRasterBand(i+1);
59 poBand->ComputeRasterMinMax(FALSE,bandMINMAX);
60 for(int j = 0;j<m_iHeigth;j++){
61 double *readBuff = new double[m_iWidth];
62 unsigned char*writeBuff = new unsigned char[m_iWidth];
63 poBand->RasterIO(GF_Read,0,j,m_iWidth,1,readBuff,m_iWidth,1,GDT_Float64,0,0);
64 int k = 0;
65 while(k<m_iWidth){
66 int tmp = 0;
67 if(eDT == GDT_Float32 || eDT == GDT_Float64){
68 tmp = int((m_grayLevel-1+0.99999)*(readBuff[k] - bandMINMAX[0])/(bandMINMAX[1] - bandMINMAX[0]));
69 }else{
70 tmp = int((m_grayLevel)*(readBuff[k] - bandMINMAX[0] - 1)/(bandMINMAX[1] - bandMINMAX[0]));
71 }
72 writeBuff[k] = unsigned char(tmp);
73 k++;
74 }
75 CPLErr str = poNewBand->RasterIO(GF_Write,0,j,m_iWidth,1,writeBuff,m_iWidth,1,GDT_Byte,0,0);
76 delete []readBuff;readBuff = NULL;
77 delete []writeBuff;writeBuff = NULL;
78 }
79 delete []bandMINMAX;bandMINMAX = NULL;
80 }
81 GDALClose((GDALDatasetH) poDS);
82 GDALClose((GDALDatasetH) poOutDS);
83 return TRUE;
84 }
85
86 double* CCo_occurrenceTextureExtration::cdf( int *h,int length )
87 {
88 int n = 0;
89 for( int i = 0; i < length - 1; i++ )
90 {
91 n += h[i];
92 }
93 double* p = new double[length];
94 int c = h[0];
95 p[0] = ( double )c / n;
96 for( int i = 1; i < length - 1; i++ )
97 {
98 c += h[i];
99 p[i] = ( double )c / n;
100 }
101 return p;
102 }
103
104 char* CCo_occurrenceTextureExtration::LoadProgSource( const char* cFilename, const char* cPreamble, size_t* szFinalLength )
105 {
106 FILE* pFileStream = NULL;
107 size_t szSourceLength;
108
109 // open the OpenCL source code file
110 pFileStream = fopen(cFilename, "rb");
111 if(pFileStream == 0)
112 {
113 return NULL;
114 }
115
116 size_t szPreambleLength = strlen(cPreamble);
117
118 // get the length of the source code
119 fseek(pFileStream, 0, SEEK_END);
120 szSourceLength = ftell(pFileStream);
121 fseek(pFileStream, 0, SEEK_SET);
122
123 // allocate a buffer for the source code string and read it in
124 char* cSourceString = (char *)malloc(szSourceLength + szPreambleLength + 1);
125 memcpy(cSourceString, cPreamble, szPreambleLength);
126 if (fread((cSourceString) + szPreambleLength, szSourceLength, 1, pFileStream) != 1)
127 {
128 fclose(pFileStream);
129 free(cSourceString);
130 return 0;
131 }
132
133 // close the file and return the total length of the combined (preamble + source) string
134 fclose(pFileStream);
135 if(szFinalLength != 0)
136 {
137 *szFinalLength = szSourceLength + szPreambleLength;
138 }
139 cSourceString[szSourceLength + szPreambleLength] = '\0';
140
141 return cSourceString;
142 }
143
144 bool CCo_occurrenceTextureExtration::GetGrayCoocurrenceMatrix(std::vector<std::vector<int>>grayValue,
145 std::vector<std::vector<float>> &m_grayMatrix,
146 int &m_count)
147 {
148 int i,j,r,c;
149 m_count = 0;
150 int ii = 0;
151
152 switch(m_AngleMode)
153 {
154 case DEGREE_0:
155 i = 0;
156 while(i<m_winHeigth){
157 j = 0;
158 while(j<m_winWidth){
159 int endR = i;
160 int endC = j + m_dis;
161 if(endC < m_winWidth){
162 r = grayValue[i][j];
163 c = grayValue[endR][endC];
164 m_grayMatrix[r][c] += 1;
165 m_grayMatrix[c][r] += 1;
166 m_count = m_count+2;
167 }
168 j++;
169 }
170 i++;
171 }
172 break;
173 case DEGREE_45:
174 i = 0;
175 while(i<m_winHeigth){
176 j = 0;
177 while(j<m_winWidth){
178 int endR = i - m_dis;
179 int endC = j + m_dis;
180 if(endR>=0 && endR < m_winHeigth && endC >=0 && endC < m_winWidth){
181 r = grayValue[i][j];
182 c = grayValue[endR][endC];
183 m_grayMatrix[r][c] += 1;
184 m_grayMatrix[c][r] += 1;
185 m_count = m_count+2;
186 }
187 j++;
188 }
189 i++;
190 }
191 break;
192 case DEGREE_90:
193 i= 0;
194 while(i<m_winHeigth){
195 j = 0;
196 while(j < m_winWidth){
197 int endR = i - m_dis;
198 int endC = j;
199 if(endR >= 0){
200 r = grayValue[i][j];
201 c = grayValue[endR][endC];
202 m_grayMatrix[r][c] += 1;
203 m_grayMatrix[c][r] += 1;
204 m_count = m_count+2;
205 }
206 j++;
207 }
208 i++;
209 }
210 break;
211 case DEGREE_135:
212 i = 0;
213 while(i<m_winHeigth){
214 j = 0;
215 while(j<m_winWidth){
216 int endR = i - m_dis;
217 int endC = j - m_dis;
218 if(endR >= 0 && endC >= 0){
219 r = grayValue[i][j];
220 c = grayValue[endR][endC];
221 m_grayMatrix[r][c] += 1;
222 m_grayMatrix[c][r] += 1;
223 m_count = m_count+2;
224 }
225 j++;
226 }
227 i++;
228 }
229 break;
230 default:
231 return FALSE;
232 }
233 return TRUE;
234 }
235
236 float CCo_occurrenceTextureExtration::CalGLCMStatistics(std::vector<std::vector<float>> m_grayMatrix,
237 const int m_count)
238 {
239 int i,j,k;
240 float imean = 0.0f;
241 float jmean = 0.0f;
242 float ivar = 0.0f;
243 float jvar = 0.0f;;
244 float res = 0.0f;
245 switch(m_TextureMode)
246 {
247 case GLCM_HOM: //同质性
248 i = 0;
249 while(i < m_grayLevel){
250 j = 0;
251 while(j < m_grayLevel){
252 if(m_grayMatrix[i][j] != 0){
253 float Pij = m_grayMatrix[i][j]*1.0f/m_count;
254 res += Pij/(1 + (i-j)*(i-j));
255 }
256 j++;
257 }
258 i++;
259 }
260 break;
261 case GLCM_CON: //对比度
262 i = 0;
263 while(i < m_grayLevel){
264 j = 0;
265 while(j < m_grayLevel){
266 if(m_grayMatrix[i][j] != 0){
267 float Pij = m_grayMatrix[i][j]*1.0f/m_count;
268 res += Pij*(i-j)*(i-j);
269 }
270 j++;
271 }
272 i++;
273 }
274 break;
275 case GLCM_DIS: //非相似性
276 i = 0;
277 while(i < m_grayLevel){
278 j = 0;
279 while(j < m_grayLevel){
280 if(m_grayMatrix[i][j] != 0){
281 float Pij = m_grayMatrix[i][j]*1.0f/m_count;
282 res += Pij*std::abs(i-j);
283 }
284 j++;
285 }
286 i++;
287 }
288 break;
289 case GLCM_COR: // 相关性
290 i = 0;
291 while(i < m_grayLevel){
292 j = 0;
293 while(j < m_grayLevel){
294 if(m_grayMatrix[i][j] != 0){
295 float Pij = m_grayMatrix[i][j]*1.0f/m_count;
296 imean += Pij * i;
297 jmean += Pij * j;
298 }
299 j++;
300 }
301 i++;
302 }
303 // var
304 i = 0;
305 while(i < m_grayLevel){
306 j = 0;
307 while(j < m_grayLevel){
308 if(m_grayMatrix[i][j] != 0){
309 float Pij = m_grayMatrix[i][j]*1.0f/m_count;
310 ivar += Pij * (i - imean)*(i - imean);
311 jvar += Pij * (j - imean)*(j - imean);
312 }
313 j++;
314 }
315 i++;
316 }
317 // xiangguanx
318 i = 0;
319 while(i < m_grayLevel){
320 j = 0;
321 while(j < m_grayLevel){
322 if(m_grayMatrix[i][j] != 0){
323 float Pij = m_grayMatrix[i][j]*1.0f/m_count;
324 res += (Pij *(i - imean)*(j - jmean))/sqrt(ivar*jvar);
325 }
326 j++;
327 }
328 i++;
329 }
330 break;
331 case GLCM_MEAN: // 均值
332 i = 0;
333 while(i < m_grayLevel){
334 j = 0;
335 while(j < m_grayLevel){
336 if(m_grayMatrix[i][j] != 0){
337 float Pij = m_grayMatrix[i][j]*1.0f/m_count;
338 imean += Pij * i;
339 //jmean += Pij * j;
340 }
341 j++;
342 }
343 i++;
344 }
345 res = imean;
346 break;
347 case GLCM_VAR: // 方差
348 // mean
349 i = 0;
350 while(i < m_grayLevel){
351 j = 0;
352 while(j < m_grayLevel){
353 if(m_grayMatrix[i][j] != 0){
354 float Pij = m_grayMatrix[i][j]*1.0f/m_count;
355 imean += Pij * i;
356 }
357 j++;
358 }
359 i++;
360 }
361 // var
362 i = 0;
363 while(i < m_grayLevel){
364 j = 0;
365 while(j < m_grayLevel){
366 if(m_grayMatrix[i][j] != 0){
367 float Pij = m_grayMatrix[i][j]*1.0f/m_count;
368 ivar += Pij * (i - imean)*(i - imean);
369 }
370 j++;
371 }
372 i++;
373 }
374 //res = (ivar + jvar)/2;
375 res = ivar;
376 break;
377 case GLCM_ASM: // 角二阶矩
378 i = 0;
379 while(i < m_grayLevel){
380 j = 0;
381 while(j < m_grayLevel){
382 if(m_grayMatrix[i][j] != 0){
383 float Pij = m_grayMatrix[i][j]*1.0f/m_count;
384 res += Pij*Pij;
385 }
386 j++;
387 }
388 i++;
389 }
390 break;
391 case GLCM_ENTROPY: // 熵
392 i = 0;
393 while(i < m_grayLevel){
394 j = 0;
395 while(j < m_grayLevel){
396 if(m_grayMatrix[i][j] != 0){
397 float Pij = m_grayMatrix[i][j]*1.0f/m_count;
398 res += -1.0f*Pij * log(Pij);
399 }
400 j++;
401 }
402 i++;
403 }
404 break;
405 case GLDV_MEAN:
406 k = 0;
407 while(k<m_grayLevel){
408 i = 0;
409 while(i < m_grayLevel){
410 j = 0;
411 while(j < m_grayLevel){
412 if(m_grayMatrix[i][j] != 0 && std::abs(i-j)== k){
413 float Pij = m_grayMatrix[i][j]*1.0f/m_count;
414 res += Pij*k;
415 }
416 j++;
417 }
418 i++;
419 }
420 k++;
421 }
422 break;
423 case GLDV_CON:
424 k = 0;
425 while(k<m_grayLevel){
426 i = 0;
427 while(i < m_grayLevel){
428 j = 0;
429 while(j < m_grayLevel){
430 if(m_grayMatrix[i][j] != 0 && std::abs(i-j)== k){
431 float Pij = m_grayMatrix[i][j]*1.0f/m_count;
432 res += Pij*k*k;
433 }
434 j++;
435 }
436 i++;
437 }
438 k++;
439 }
440 break;
441 case GLDV_ASM:
442 k = 0;
443 while(k<m_grayLevel){
444 i = 0;
445 while(i < m_grayLevel){
446 j = 0;
447 while(j < m_grayLevel){
448 if(m_grayMatrix[i][j] != 0 && std::abs(i-j)== k){
449 float Pij = m_grayMatrix[i][j]*1.0f/m_count;
450 res += Pij*Pij;
451 }
452 j++;
453 }
454 i++;
455 }
456 k++;
457 }
458 break;
459 }
460 return res;
461 }
462
463 bool CCo_occurrenceTextureExtration::ExtraTextureWinEXE()
464 {
465 //读影像
466 GDALAllRegister();
467 GDALDataset *poDS = (GDALDataset*)GDALOpen(m_strReduceLevelFileName.c_str(),GA_ReadOnly);
468 if(!poDS){
469 return FALSE;
470 }
471 //GDALRasterBand *poBand = poDS->GetRasterBand(1);
472 double adfMSSGeoTransform[6] = {0};
473 poDS->GetGeoTransform(adfMSSGeoTransform);
474 //创建文件
475 GDALDriver *poOutDriver = (GDALDriver*)GDALGetDriverByName("GTIFF");
476 GDALDataset *poOutDS = poOutDriver->Create(m_strOutFileName.c_str(),
477 m_iWidth,m_iHeigth,m_iBandCount,GDT_Float32,NULL);
478 poOutDS->SetGeoTransform(adfMSSGeoTransform);
479 poOutDS->SetProjection(poDS->GetProjectionRef());
480 int b = 0;
481 while(b<m_iBandCount){
482 GDALRasterBand *poBand = poDS->GetRasterBand(b+1);
483 GDALRasterBand *poNewBand = poOutDS->GetRasterBand(b+1);
484 int i = 0;
485 while(i<m_iHeigth){
486 int startR = i;
487 int endR = startR + m_winHeigth - 1;
488 if(endR > m_iHeigth-1){
489 endR = m_iHeigth - 1;
490 startR = endR - m_winHeigth + 1;
491 }
492 unsigned char *readBuff = new unsigned char[m_iWidth * m_winHeigth];
493 float *writeBuff = new float[m_iWidth* m_winHeigth];
494 poBand->RasterIO(GF_Read,0,startR,m_iWidth,m_winHeigth,readBuff,m_iWidth,
495 m_winHeigth,GDT_Byte,0,0);
496
497 int j = 0;
498 while(j < m_iWidth){
499 int startC = j;
500 int endC = startC + m_winWidth - 1;
501 if(endC > m_iWidth - 1){
502 endC = m_iWidth - 1;
503 startC = endC - m_winWidth + 1;
504 }
505 std::vector<std::vector<int>> grayValue;
506 grayValue = std::vector<std::vector<int>>(m_winHeigth,std::vector<int>(m_winWidth,0));
507 std::vector<std::vector<float>> grayMatrix;
508 grayMatrix = std::vector<std::vector<float>>(m_grayLevel,std::vector<float>(m_grayLevel,0.0f));
509 int k = 0;
510 while(k < m_winHeigth){
511 int l = 0 ;
512 while(l < m_winWidth){
513 grayValue[k][l] = int(readBuff[k*m_iWidth +startC+l]);
514 l++;
515 }
516 k++;
517 }
518 int count = 0;
519 if(!GetGrayCoocurrenceMatrix(grayValue,grayMatrix,count)){
520 return FALSE;
521 }
522 float tmp = CalGLCMStatistics(grayMatrix,count);
523 k = 0;
524 while(k < m_winHeigth){
525 int l = 0 ;
526 while(l < m_winWidth){
527 writeBuff[k*m_iWidth +startC+l] = tmp;
528 l++;
529 }
530 k++;
531 }
532 j = j + m_winWidth;
533 }
534 poNewBand->RasterIO(GF_Write,0,startR,m_iWidth,m_winHeigth,writeBuff,m_iWidth,
535 m_winHeigth,GDT_Float32,0,0);
536 //更新灰度共生矩阵
537 delete []readBuff;readBuff = NULL;
538 delete []writeBuff;writeBuff = NULL;
539 i = i + m_winHeigth;
540 }
541 b++;
542 }
543 GDALClose((GDALDatasetH)poDS);
544 GDALClose((GDALDatasetH)poOutDS);
545 return TRUE;
546 }
547
548 bool CCo_occurrenceTextureExtration::ExtraTexturePixEXE()
549 {
550 //读影像
551 GDALAllRegister();
552 GDALDataset *poDS = (GDALDataset*)GDALOpen(m_strReduceLevelFileName.c_str(),GA_ReadOnly);
553 if(!poDS){
554 return FALSE;
555 }
556 double adfMSSGeoTransform[6] = {0};
557 poDS->GetGeoTransform(adfMSSGeoTransform);
558 //创建文件
559 GDALDriver *poOutDriver = (GDALDriver*)GDALGetDriverByName("GTIFF");
560 GDALDataset *poOutDS = poOutDriver->Create(m_strOutFileName.c_str(),
561 m_iWidth,m_iHeigth,m_iBandCount,GDT_Float32,NULL);
562 poOutDS->SetGeoTransform(adfMSSGeoTransform);
563 poOutDS->SetProjection(poDS->GetProjectionRef());
564 int nLineSpace,nPixSpace,nBandSpace;
565 nLineSpace = sizeof(unsigned char)*m_iWidth*m_iBandCount;
566 nPixSpace = 0;
567 nBandSpace = sizeof(unsigned char)*m_iWidth;
568
569 int iNumRow = 128; // 分块读取的行数
570 int overloadPix = m_winHeigth-1;
571 if(m_iHeigth < iNumRow && m_iHeigth > m_winHeigth){
572 iNumRow = m_winHeigth;
573 }
574 if(m_iHeigth < m_winHeigth){
575 return FALSE;
576 }
577 int loopNum = (m_iHeigth + iNumRow - 1)/(iNumRow - overloadPix);
578 float stepProcess = 100.0f/(100.0f*loopNum*(iNumRow - overloadPix));
579
580 for(int i = 0;i<loopNum;i++){
581 int tmpRowNum = iNumRow;
582 int startR = i*iNumRow - i*(overloadPix/2);
583 if(startR < 0)
584 {
585 startR = 0;
586 }
587 int endR = startR + tmpRowNum - 1;
588 if(endR > m_iHeigth -1 ){
589 endR = m_iHeigth -1;
590 tmpRowNum = endR - startR + 1;
591 }
592 unsigned char *readBuff = new unsigned char[tmpRowNum*m_iWidth*m_iBandCount];
593 float *writeBuff = new float[m_iWidth*m_iBandCount*(tmpRowNum - overloadPix)];
594 poDS->RasterIO(GF_Read,0,startR,m_iWidth,tmpRowNum,readBuff,m_iWidth,
595 tmpRowNum,GDT_Byte,m_iBandCount,m_BandMap,nPixSpace,nLineSpace,nBandSpace);
596
597 for(int iR = overloadPix/2;iR < tmpRowNum - overloadPix/2;iR++)
598 {
599 // zhuyong 2016 10 17
600 #pragma omp parallel for num_threads(2)
601 for(int j = 0;j < m_iWidth;j++)
602 {
603 int startC = j - m_winWidth/2;
604 int endC = j + m_winWidth/2;
605 if(startC<0){
606 startC = 0;
607 endC = startC + m_winWidth -1;
608 }
609 if(endC > m_iWidth-1){
610 endC = m_iWidth -1;
611 startC = endC - m_winWidth +1;
612 }
613 // 波段计算
614 int b = 0;
615 while(b < m_iBandCount){
616 // 获取灰度值
617 std::vector<std::vector<int>> grayValue;
618 grayValue = std::vector<std::vector<int>>(m_winHeigth,std::vector<int>(m_winWidth,0));
619 std::vector<std::vector<float>> grayMatrix;
620 grayMatrix = std::vector<std::vector<float>>(m_grayLevel,std::vector<float>(m_grayLevel,0.0f));
621 int k = 0;
622 while(k < m_winHeigth){
623 int l = 0 ;
624 while(l < m_winWidth){
625 grayValue[k][l] = int(readBuff[(iR - overloadPix/2 + k)*m_iWidth*m_iBandCount + b*m_iWidth + startC+l]);
626 l++;
627 }
628 k++;
629 }
630 int count = 0;
631 if(!GetGrayCoocurrenceMatrix(grayValue,grayMatrix,count)){
632 }
633 float tmp = CalGLCMStatistics(grayMatrix,count);
634 writeBuff[(iR- overloadPix/2)*m_iWidth*m_iBandCount+ b*m_iWidth + j] = tmp;
635 b++;
636 }
637 }
638 }
639 poOutDS->RasterIO(GF_Write,0,startR+overloadPix/2,m_iWidth,tmpRowNum - overloadPix,writeBuff,
640 m_iWidth,tmpRowNum - overloadPix,GDT_Float32,m_iBandCount,m_BandMap,nPixSpace,sizeof(float)*m_iWidth*m_iBandCount,
641 sizeof(float)*m_iWidth);
642 delete []readBuff;readBuff = NULL;
643 delete []writeBuff;writeBuff = NULL;
644 }
645
646 GDALClose((GDALDatasetH)poDS);
647 GDALClose((GDALDatasetH)poOutDS);
648 return TRUE;
649 }
650
651 int CCo_occurrenceTextureExtration::ExtraTexturePixOpenCLEXE()
652 {
653 //读影像
654 GDALAllRegister();
655 GDALDataset *poDS = (GDALDataset*)GDALOpen(m_strReduceLevelFileName.c_str(),GA_ReadOnly);
656 if(!poDS){
657 return FALSE;
658 }
659 double adfMSSGeoTransform[6] = {0};
660 poDS->GetGeoTransform(adfMSSGeoTransform);
661 //创建文件
662 GDALDriver *poOutDriver = (GDALDriver*)GDALGetDriverByName("GTIFF");
663 GDALDataset *poOutDS = poOutDriver->Create(m_strOutFileName.c_str(),
664 m_iWidth,m_iHeigth,m_iBandCount,GDT_Float32,NULL);
665 poOutDS->SetGeoTransform(adfMSSGeoTransform);
666 poOutDS->SetProjection(poDS->GetProjectionRef());
667
668 int nLineSpace,nPixSpace,nBandSpace;
669 nLineSpace = sizeof(unsigned char)*m_iWidth*m_iBandCount;
670 nPixSpace = 0;
671 nBandSpace = sizeof(unsigned char)*m_iWidth;
672
673 int iNumRow = 32; // 分块读取的行数
674 int overloadPix = m_winHeigth-1;
675 if(m_iHeigth < iNumRow && m_iHeigth > m_winHeigth){
676 iNumRow = m_winHeigth;
677 }
678 if(m_iHeigth < m_winHeigth){
679 return FALSE;
680 }
681 int loopNum = ceil((m_iHeigth*1.0f - iNumRow)/(iNumRow - overloadPix)) + 1;
682
683 // OpenCL部分 =============== 1 创建平台
684 cl_uint num_platforms;
685 cl_int ret = clGetPlatformIDs(0,NULL,&num_platforms);
686 if(ret != CL_SUCCESS || num_platforms < 1){
687 printf("clGetPlatformIDs Error\n");
688 return -1;
689 }
690 cl_platform_id platform_id = NULL;
691 ret = clGetPlatformIDs(1,&platform_id,NULL);
692 if(ret != CL_SUCCESS){
693 printf("clGetPlatformIDs Error2\n");
694 return -1;
695 }
696
697 // OpenCL部分 =============== 2 获得设备
698 cl_uint num_devices;
699 ret = clGetDeviceIDs(platform_id,CL_DEVICE_TYPE_GPU,0,NULL,
700 &num_devices);
701 if(ret != CL_SUCCESS || num_devices < 1){
702 printf("clGetDeviceIDs Error\n");
703 return -1;
704 }
705 cl_device_id device_id;
706 ret = clGetDeviceIDs(platform_id,CL_DEVICE_TYPE_GPU,1,&device_id,NULL);
707 if(ret != CL_SUCCESS){
708 printf("clGetDeviceIDs Error2\n");
709 return -1;
710 }
711
712 // OpenCL部分 =============== 3 创建Context
713 cl_context_properties props[] = {CL_CONTEXT_PLATFORM,
714 (cl_context_properties)platform_id,0};
715 cl_context context = NULL;
716 context = clCreateContext(props,1,&device_id,NULL,NULL,&ret);
717 if(ret != CL_SUCCESS || context == NULL){
718 printf("clCreateContext Error\n");
719 return -1;
720 }
721
722 // OpenCL部分 =============== 4 创建Command Queue
723 cl_command_queue command_queue = NULL;
724 command_queue = clCreateCommandQueue(context,device_id,0,&ret);
725 if(ret != CL_SUCCESS || command_queue == NULL){
726 printf("clCreateCommandQueue Error\n");
727 return -1;
728 }
729
730 // OpenCL部分 =============== 6 创建编译Program
731 const char *strfile = "D:\\PIE3\\src\\Test\\TextOpecCLResample\\TextOpecCLResample\\Co_occurrenceMatrixKernel.txt";
732 size_t lenSource = 0;
733 char *kernelSource = LoadProgSource(strfile,"",&lenSource);
734 cl_program *programs = (cl_program *)malloc(loopNum*sizeof(cl_program));
735 memset(programs,0,sizeof(cl_program)*loopNum);
736
737 cl_kernel *kernels = (cl_kernel*)malloc(loopNum*sizeof(cl_kernel));
738 memset(kernels,0,sizeof(cl_kernel)*loopNum);
739
740 int startR = 0;
741 for(int i = 0;i<loopNum;i++)
742 {
743 int tmpRowNum = iNumRow;
744 if(startR < 0)
745 {
746 startR = 0;
747 }
748 int endR = startR + tmpRowNum - 1;
749 if(endR > m_iHeigth -1 ){
750 endR = m_iHeigth -1;
751 tmpRowNum = endR - startR + 1;
752 }
753 unsigned char *readBuff = new unsigned char[tmpRowNum*m_iWidth*m_iBandCount];
754 float *writeBuff = new float[m_iWidth*m_iBandCount*(tmpRowNum - overloadPix)];
755 poDS->RasterIO(GF_Read,0,startR,m_iWidth,tmpRowNum,readBuff,m_iWidth,
756 tmpRowNum,GDT_Byte,m_iBandCount,m_BandMap,nPixSpace,nLineSpace,nBandSpace);
757
758 // OpenCL部分 =============== 5 创建Memory Object
759 cl_mem mem_read = NULL;
760 mem_read = clCreateBuffer(context,CL_MEM_READ_WRITE | CL_MEM_USE_HOST_PTR,
761 sizeof(cl_uchar)*tmpRowNum*m_iWidth*m_iBandCount,readBuff,&ret);
762 if(ret != CL_SUCCESS || NULL == mem_read){
763 printf("clCreateBuffer Error\n");
764 return -1;
765 }
766
767 cl_mem mem_write = NULL;
768 mem_write = clCreateBuffer(context,CL_MEM_READ_WRITE | CL_MEM_USE_HOST_PTR,
769 sizeof(cl_float)*m_iWidth*m_iBandCount*(tmpRowNum - overloadPix),writeBuff,&ret);
770 if(ret != CL_SUCCESS || NULL == mem_write){
771 printf("clCreateBuffer Error\n");
772 return -1;
773 }
774
775 // OpenCL部分 =============== 6 创建编译Program
776 //const char *strfile = "D:\\PIE3\\src\\Test\\TextOpecCLResample\\TextOpecCLResample\\Co_occurrenceMatrixKernel.txt";
777 //size_t lenSource = 0;
778 //char *kernelSource = LoadProgSource(strfile,"",&lenSource);
779 //cl_program program = NULL;
780 programs[i] = clCreateProgramWithSource(context,1,(const char**)&kernelSource,
781 NULL,&ret);
782 if(ret != CL_SUCCESS || NULL == programs[i]){
783 printf("clCreateProgramWithSource Error\n");
784 return -1;
785 }
786 ret = clBuildProgram(programs[i],1,&device_id,NULL,NULL,NULL);
787 if(ret != CL_SUCCESS){
788 char* build_log;
789 size_t log_size;
790 //查询日志的大小
791 clGetProgramBuildInfo(programs[i], device_id, CL_PROGRAM_BUILD_LOG, 0, NULL, &log_size);
792 build_log = new char[log_size+1];
793 //获得编译日志信息
794 ret = clGetProgramBuildInfo(programs[i], device_id, CL_PROGRAM_BUILD_LOG, log_size, build_log, NULL);
795 build_log[log_size] = '\0';
796 printf("%s\n",build_log);
797 printf("编译失败!");
798 delete []build_log;
799 return -1;
800 }
801
802 // OpenCL部分 =============== 7 创建Kernel
803 //cl_kernel kernel = NULL;
804 kernels[i] = clCreateKernel(programs[i],"Co_occurrenceMatrixKernel",&ret);
805 if(ret != CL_SUCCESS || NULL == kernels[i]){
806 printf("clCreateProgramWithSource Error\n");
807 return -1;
808 }
809
810 // OpenCL部分 =============== 8 设置Kernel参数
811 ret = clSetKernelArg(kernels[i],0,sizeof(cl_mem),&mem_read);
812 ret |= clSetKernelArg(kernels[i],1,sizeof(cl_mem),&mem_write);
813 ret |= clSetKernelArg(kernels[i],2,sizeof(cl_int),&m_iHeigth);
814 ret |= clSetKernelArg(kernels[i],3,sizeof(cl_int),&m_iWidth);
815 ret |= clSetKernelArg(kernels[i],4,sizeof(cl_int),&m_iBandCount);
816 ret |= clSetKernelArg(kernels[i],5,sizeof(cl_int),&tmpRowNum);
817 ret |= clSetKernelArg(kernels[i],6,sizeof(cl_int),&overloadPix);
818 ret |= clSetKernelArg(kernels[i],7,sizeof(cl_int),&m_grayLevel);
819 ret |= clSetKernelArg(kernels[i],8,sizeof(cl_int),&m_winHeigth);
820 ret |= clSetKernelArg(kernels[i],9,sizeof(cl_int),&m_winWidth);
821 ret |= clSetKernelArg(kernels[i],10,sizeof(cl_int),&m_dis);
822 ret |= clSetKernelArg(kernels[i],11,sizeof(cl_int),&m_AngleMode);
823 ret |= clSetKernelArg(kernels[i],12,sizeof(cl_int),&m_TextureMode);
824 if(ret != CL_SUCCESS){
825 printf("clSetKernelArg Error\n");
826 return -1;
827 }
828
829 // OpenCL部分 =============== 9 设置Group Size
830 cl_uint work_dim = 2;
831 size_t global_work_size[] = {m_iWidth,tmpRowNum};
832 size_t *local_work_size = NULL;
833
834 // OpenCL部分 =============== 10 执行内核
835 ret = clEnqueueNDRangeKernel(command_queue,kernels[i],work_dim,NULL,global_work_size,
836 local_work_size,0,NULL,NULL);
837 ret |= clFinish(command_queue);
838 if(ret != CL_SUCCESS){
839 printf("clEnqueueNDRangeKernel Error\n");
840 return -1;
841 }
842
843 writeBuff = (float*)clEnqueueMapBuffer(command_queue,mem_write,CL_TRUE,CL_MAP_READ | CL_MAP_WRITE,
844 0,sizeof(cl_float)*m_iWidth*m_iBandCount*(tmpRowNum - overloadPix),0,NULL,NULL,&ret);
845 //ret = clEnqueueReadBuffer(command_queue,mem_resample,CL_TRUE,0,
846 // sizeof(cl_float)*tmpRowNum*resampleWidth*mssBandCount,(void*)resampleBuf,0,NULL,NULL);
847 if(ret != CL_SUCCESS){
848 printf("clEnqueueMapBuffer Error\n");
849 return -1;
850 }
851 poOutDS->RasterIO(GF_Write,0,startR+overloadPix/2,m_iWidth,tmpRowNum - overloadPix,writeBuff,
852 m_iWidth,tmpRowNum - overloadPix,GDT_Float32,m_iBandCount,m_BandMap,nPixSpace,sizeof(float)*m_iWidth*m_iBandCount,
853 sizeof(float)*m_iWidth);
854 delete []readBuff;readBuff = NULL;
855 delete []writeBuff;writeBuff = NULL;
856
857 // OpenCL部分 =============== 12 释放资源
858 if(NULL != mem_read) clReleaseMemObject(mem_read);
859 if(NULL != mem_write) clReleaseMemObject(mem_write);
860 std::cout<<i<<" "<<loopNum<<std::endl;
861 startR = startR + tmpRowNum - overloadPix;
862 }
863 // OpenCL部分 =============== 12 释放资源
864 int i = 0;
865 while(i < loopNum){
866 if(NULL != kernels[i]) clReleaseKernel(kernels[i]);
867 if(NULL != programs[i]) clReleaseProgram(programs[i]);
868 i++;
869 }
870 if(NULL != command_queue) clReleaseCommandQueue(command_queue);
871 if(NULL != context) clReleaseContext(context);
872 GDALClose((GDALDatasetH)poDS);
873 GDALClose((GDALDatasetH)poOutDS);
874 return 0;
875 }GPU核函数：
1 #pragma OPENCL EXTENSION cl_amd_printf:enable
2
3 int GetGrayCoocurrenceMatrix(const int *grayValue,
4 float *grayMatrix,
5 int winHeight,
6 int winWidth,
7 int grayLevel,
8 int dis,
9 int AngleMode)
10 {
11 int i,j,r,c;
12 int count = 0;
13 switch(AngleMode)
14 {
15 case 0:
16 i = 0;
17 while(i<winHeight){
18 j = 0;
19 while(j<winWidth){
20 int endR = i;
21 int endC = j + dis;
22 if(endC < winWidth){
23 r = grayValue[i*winHeight + j];
24 c = grayValue[endR*winHeight + endC];
25 grayMatrix[r*grayLevel + c] += 1;
26 grayMatrix[c*grayLevel + r] += 1;
27 count = count+2;
28 }
29 j++;
30 }
31 i++;
32 }
33 break;
34 case 1:
35 i = 0;
36 while(i<winHeight){
37 j = 0;
38 while(j<winWidth){
39 int endR = i - dis;
40 int endC = j + dis;
41 if(endR>=0 && endR < winHeight && endC >=0 && endC < winWidth){
42 r = grayValue[i*winHeight + j];
43 c = grayValue[endR*winHeight + endC];
44 grayMatrix[r*grayLevel + c] += 1;
45 grayMatrix[c*grayLevel + r] += 1;
46 count = count+2;
47 }
48 j++;
49 }
50 i++;
51 }
52 break;
53 case 2:
54 i= 0;
55 while(i<winHeight){
56 j = 0;
57 while(j < winWidth){
58 int endR = i - dis;
59 int endC = j;
60 if(endR >= 0){
61 r = grayValue[i*winHeight + j];
62 c = grayValue[endR*winHeight + endC];
63 grayMatrix[r*grayLevel + c] += 1;
64 grayMatrix[c*grayLevel + r] += 1;
65 count = count+2;
66 }
67 j++;
68 }
69 i++;
70 }
71 break;
72 case 3:
73 i = 0;
74 while(i<winHeight){
75 j = 0;
76 while(j<winWidth){
77 int endR = i - dis;
78 int endC = j - dis;
79 if(endR >= 0 && endC >= 0){
80 r = grayValue[i*winHeight + j];
81 c = grayValue[endR*winHeight + endC];
82 grayMatrix[r*grayLevel + c] += 1;
83 grayMatrix[c*grayLevel + r] += 1;
84 count = count+2;
85 }
86 j++;
87 }
88 i++;
89 }
90 break;
91 default:
92 return 0;
93 }
94 return count;
95 }
96
97 float CalGLCMStatistics(float *grayMatrix,
98 const int grayLevel,
99 const int count,
100 int TextureMode)
101 {
102 int i,j,k;
103 float imean = 0.0f;
104 float jmean = 0.0f;
105 float ivar = 0.0f;
106 float jvar = 0.0f;;
107 float res = 0.0f;
108 switch(TextureMode)
109 {
110 case 0: //同质性
111 i = 0;
112 while(i < grayLevel){
113 j = 0;
114 while(j < grayLevel){
115 if(grayMatrix[i*grayLevel + j] != 0){
116 float Pij = grayMatrix[i*grayLevel + j]*1.0f/count;
117 res += Pij/(1 + (i-j)*(i-j));
118 }
119 j++;
120 }
121 i++;
122 }
123 break;
124 case 1: //对比度
125 i = 0;
126 while(i < grayLevel){
127 j = 0;
128 while(j < grayLevel){
129 if(grayMatrix[i*grayLevel + j] != 0){
130 float Pij = grayMatrix[i*grayLevel + j]*1.0f/count;
131 res += Pij*(i-j)*(i-j);
132 }
133 j++;
134 }
135 i++;
136 }
137 break;
138 case 2: //非相似性
139 i = 0;
140 while(i < grayLevel){
141 j = 0;
142 while(j < grayLevel){
143 if(grayMatrix[i*grayLevel + j] != 0){
144 float Pij = grayMatrix[i*grayLevel + j]*1.0f/count;
145 res += Pij*abs(i-j);
146 }
147 j++;
148 }
149 i++;
150 }
151 break;
152 case 3: // 相关性
153 i = 0;
154 while(i < grayLevel){
155 j = 0;
156 while(j < grayLevel){
157 if(grayMatrix[i*grayLevel + j] != 0){
158 float Pij = grayMatrix[i*grayLevel + j]*1.0f/count;
159 imean += Pij * i;
160 jmean += Pij * j;
161 }
162 j++;
163 }
164 i++;
165 }
166 // var
167 i = 0;
168 while(i < grayLevel){
169 j = 0;
170 while(j < grayLevel){
171 if(grayMatrix[i*grayLevel + j] != 0){
172 float Pij = grayMatrix[i*grayLevel + j]*1.0f/count;
173 ivar += Pij * (i - imean)*(i - imean);
174 jvar += Pij * (j - imean)*(j - imean);
175 }
176 j++;
177 }
178 i++;
179 }
180 // xiangguanx
181 i = 0;
182 while(i < grayLevel){
183 j = 0;
184 while(j < grayLevel){
185 if(grayMatrix[i*grayLevel + j] != 0){
186 float Pij = grayMatrix[i*grayLevel + j]*1.0f/count;
187 res += (Pij *(i - imean)*(j - jmean))/sqrt(ivar*jvar);
188 }
189 j++;
190 }
191 i++;
192 }
193 break;
194 case 4: // 均值
195 i = 0;
196 while(i < grayLevel){
197 j = 0;
198 while(j < grayLevel){
199 if(grayMatrix[i*grayLevel + j] != 0){
200 float Pij = grayMatrix[i*grayLevel + j]*1.0f/count;
201 imean += Pij * i;
202 //jmean += Pij * j;
203 }
204 j++;
205 }
206 i++;
207 }
208 res = imean;
209 break;
210 case 5: // 方差
211 // mean
212 i = 0;
213 while(i < grayLevel){
214 j = 0;
215 while(j < grayLevel){
216 if(grayMatrix[i*grayLevel + j] != 0){
217 float Pij = grayMatrix[i*grayLevel + j]*1.0f/count;
218 imean += Pij * i;
219 }
220 j++;
221 }
222 i++;
223 }
224 // var
225 i = 0;
226 while(i < grayLevel){
227 j = 0;
228 while(j < grayLevel){
229 if(grayMatrix[i*grayLevel + j] != 0){
230 float Pij = grayMatrix[i*grayLevel + j]*1.0f/count;
231 ivar += Pij * (i - imean)*(i - imean);
232 }
233 j++;
234 }
235 i++;
236 }
237 res = ivar;
238 break;
239 case 6: // 角二阶矩
240 i = 0;
241 while(i < grayLevel){
242 j = 0;
243 while(j < grayLevel){
244 if(grayMatrix[i*grayLevel + j] != 0){
245 float Pij = grayMatrix[i*grayLevel + j]*1.0f/count;
246 res += Pij*Pij;
247 }
248 j++;
249 }
250 i++;
251 }
252 break;
253 case 7: // 熵
254 i = 0;
255 while(i < grayLevel){
256 j = 0;
257 while(j < grayLevel){
258 if(grayMatrix[i*grayLevel + j] != 0){
259 float Pij = grayMatrix[i*grayLevel + j]*1.0f/count;
260 res += -1.0f*Pij * log(Pij);
261 }
262 j++;
263 }
264 i++;
265 }
266 break;
267 case 8:
268 k = 0;
269 while(k<grayLevel){
270 i = 0;
271 while(i < grayLevel){
272 j = 0;
273 while(j < grayLevel){
274 if(grayMatrix[i*grayLevel + j] != 0 && abs(i-j)== k){
275 float Pij = grayMatrix[i*grayLevel + j]*1.0f/count;
276 res += Pij*k;
277 }
278 j++;
279 }
280 i++;
281 }
282 k++;
283 }
284 break;
285 case 9:
286 k = 0;
287 while(k<grayLevel){
288 i = 0;
289 while(i < grayLevel){
290 j = 0;
291 while(j < grayLevel){
292 if(grayMatrix[i*grayLevel + j] != 0 && abs(i-j)== k){
293 float Pij = grayMatrix[i*grayLevel + j]*1.0f/count;
294 res += Pij*k*k;
295 }
296 j++;
297 }
298 i++;
299 }
300 k++;
301 }
302 break;
303 case 10:
304 k = 0;
305 while(k<grayLevel){
306 i = 0;
307 while(i < grayLevel){
308 j = 0;
309 while(j < grayLevel){
310 if(grayMatrix[i*grayLevel + j] != 0 && abs(i-j)== k){
311 float Pij = grayMatrix[i*grayLevel + j]*1.0f/count;
312 res += Pij*Pij;
313 }
314 j++;
315 }
316 i++;
317 }
318 k++;
319 }
320 break;
321 }
322 return res;
323 }
324
325
326
327
328
329
330 __kernel void Co_occurrenceMatrixKernel(__global const unsigned char *oriData,
331 __global float *writeData,
332 int iHeight,
333 int iWidth,
334 int iBandCount,
335 int curRowNum,
336 int overloadPix,
337 int grayLevel,
338 int winHeight,
339 int winWidth,
340 int dis,
341 int AngleMode,
342 int TextureMode)
343 {
344 int idx = get_global_id(1); // 采样行
345 int idy = get_global_id(0); // 采样列
346 if(idx >= overloadPix/2 && idx < curRowNum - overloadPix/2 && idy < iWidth)
347 {
348 // 获得灰度值
349 int startC = idy - winWidth/2;
350 int endC = idy + winWidth/2;
351 if(startC<0){
352 startC = 0;
353 endC = startC + winWidth -1;
354 }
355 if(endC > iWidth-1){
356 endC = iWidth -1;
357 startC = endC - winWidth +1;
358 }
359 int b = 0;
360 while(b < iBandCount)
361 {
362 int grayValue[625] = {0};
363 int k = 0;
364 while(k < winHeight){
365 int l = 0 ;
366 while(l < winWidth){
367 int tmpP = (int)(idx - overloadPix/2 + k)*iWidth*iBandCount + b*iWidth + startC+l;
368 grayValue[k*winHeight + l] = (int)oriData[tmpP];
369 l++;
370 }
371 k++;
372 }
373 float grayMatrix[1024] = {0};
374 int count = GetGrayCoocurrenceMatrix(grayValue,grayMatrix,winHeight,winWidth,grayLevel,dis,AngleMode);
375 float tmp = CalGLCMStatistics(grayMatrix,grayLevel,count,TextureMode);
376 writeData[(idx- overloadPix/2)*iWidth*iBandCount+ b*iWidth + idy] = tmp;
377 b++;
378 }
379 }
380 }本文章为转载内容，我们尊重原作者对文章享有的著作权。如有内容错误或侵权问题，欢迎原作者联系我们进行内容更正或删除文章。
