diff --git a/BaseCommonLibrary/BaseTool/ImageOperatorBase.h b/BaseCommonLibrary/BaseTool/ImageOperatorBase.h index b16127d..3af8396 100644 --- a/BaseCommonLibrary/BaseTool/ImageOperatorBase.h +++ b/BaseCommonLibrary/BaseTool/ImageOperatorBase.h @@ -357,17 +357,20 @@ inline std::shared_ptr readDataArr(gdalImage& imgds, long start_row, long sta rows_count = start_row + rows_count <= imgds.height ? rows_count : imgds.height - start_row; cols_count = start_col + cols_count <= imgds.width ? cols_count : imgds.width - start_col; + int64_t pixel_count64 = static_cast(rows_count)* static_cast(cols_count); + + //Eigen::MatrixXd datamatrix(rows_count, cols_count); if (gdal_datatype == GDT_Byte) { - char* temp = new char[rows_count * cols_count]; + char* temp = new char[pixel_count64]; demBand->RasterIO(GF_Read, start_col, start_row, cols_count, rows_count, temp, cols_count, rows_count, gdal_datatype, 0, 0); - result = std::shared_ptr(new T[rows_count * cols_count], delArrPtr); + result = std::shared_ptr(new T[pixel_count64], delArrPtr); if (method == GDALREADARRCOPYMETHOD::MEMCPYMETHOD) { - std::memcpy(result.get(), temp, rows_count * cols_count); + std::memcpy(result.get(), temp, pixel_count64); } else if (method == GDALREADARRCOPYMETHOD::VARIABLEMETHOD) { - long count = rows_count * cols_count; + long count = pixel_count64; for (long i = 0; i < count; i++) { result.get()[i] = T(temp[i]); } @@ -376,14 +379,14 @@ inline std::shared_ptr readDataArr(gdalImage& imgds, long start_row, long sta delete[] temp; } else if (gdal_datatype == GDT_UInt16) { - unsigned short* temp = new unsigned short[rows_count * cols_count]; + unsigned short* temp = new unsigned short[pixel_count64]; demBand->RasterIO(GF_Read, start_col, start_row, cols_count, rows_count, temp, cols_count, rows_count, gdal_datatype, 0, 0); - result = std::shared_ptr(new T[rows_count * cols_count], delArrPtr); + result = std::shared_ptr(new T[pixel_count64], delArrPtr); if (method == GDALREADARRCOPYMETHOD::MEMCPYMETHOD) { - std::memcpy(result.get(), temp, rows_count * cols_count); + std::memcpy(result.get(), temp, pixel_count64); } else if (method == GDALREADARRCOPYMETHOD::VARIABLEMETHOD) { - long count = rows_count * cols_count; + long count = pixel_count64; for (long i = 0; i < count; i++) { result.get()[i] = T(temp[i]); } @@ -391,14 +394,14 @@ inline std::shared_ptr readDataArr(gdalImage& imgds, long start_row, long sta delete[] temp; } else if (gdal_datatype == GDT_Int16) { - short* temp = new short[rows_count * cols_count]; + short* temp = new short[pixel_count64]; demBand->RasterIO(GF_Read, start_col, start_row, cols_count, rows_count, temp, cols_count, rows_count, gdal_datatype, 0, 0); - result = std::shared_ptr(new T[rows_count * cols_count], delArrPtr); + result = std::shared_ptr(new T[pixel_count64], delArrPtr); if (method == GDALREADARRCOPYMETHOD::MEMCPYMETHOD) { - std::memcpy(result.get(), temp, rows_count * cols_count); + std::memcpy(result.get(), temp, pixel_count64); } else if (method == GDALREADARRCOPYMETHOD::VARIABLEMETHOD) { - long count = rows_count * cols_count; + long count = pixel_count64; for (long i = 0; i < count; i++) { result.get()[i] = T(temp[i]); } @@ -406,14 +409,14 @@ inline std::shared_ptr readDataArr(gdalImage& imgds, long start_row, long sta delete[] temp; } else if (gdal_datatype == GDT_UInt32) { - unsigned int* temp = new unsigned int[rows_count * cols_count]; + unsigned int* temp = new unsigned int[pixel_count64]; demBand->RasterIO(GF_Read, start_col, start_row, cols_count, rows_count, temp, cols_count, rows_count, gdal_datatype, 0, 0); - result = std::shared_ptr(new T[rows_count * cols_count], delArrPtr); + result = std::shared_ptr(new T[pixel_count64], delArrPtr); if (method == GDALREADARRCOPYMETHOD::MEMCPYMETHOD) { - std::memcpy(result.get(), temp, rows_count * cols_count); + std::memcpy(result.get(), temp, pixel_count64); } else if (method == GDALREADARRCOPYMETHOD::VARIABLEMETHOD) { - long count = rows_count * cols_count; + long count = pixel_count64; for (long i = 0; i < count; i++) { result.get()[i] = T(temp[i]); } @@ -421,14 +424,14 @@ inline std::shared_ptr readDataArr(gdalImage& imgds, long start_row, long sta delete[] temp; } else if (gdal_datatype == GDT_Int32) { - int* temp = new int[rows_count * cols_count]; + int* temp = new int[pixel_count64]; demBand->RasterIO(GF_Read, start_col, start_row, cols_count, rows_count, temp, cols_count, rows_count, gdal_datatype, 0, 0); - result = std::shared_ptr(new T[rows_count * cols_count], delArrPtr); + result = std::shared_ptr(new T[pixel_count64], delArrPtr); if (method == GDALREADARRCOPYMETHOD::MEMCPYMETHOD) { - std::memcpy(result.get(), temp, rows_count * cols_count); + std::memcpy(result.get(), temp, pixel_count64); } else if (method == GDALREADARRCOPYMETHOD::VARIABLEMETHOD) { - long count = rows_count * cols_count; + long count = pixel_count64; for (long i = 0; i < count; i++) { result.get()[i] = T(temp[i]); } @@ -436,7 +439,7 @@ inline std::shared_ptr readDataArr(gdalImage& imgds, long start_row, long sta delete[] temp; } // else if (gdal_datatype == GDT_UInt64) { - // unsigned long* temp = new unsigned long[rows_count * cols_count]; + // unsigned long* temp = new unsigned long[pixel_count64]; // demBand->RasterIO(GF_Read, start_col, start_row, cols_count, rows_count, temp, cols_count, //rows_count, gdal_datatype, 0, 0); for (int i = 0; i < rows_count; i++) { for (int j = 0; j < //cols_count; j++) { datamatrix(i, j) = temp[i * cols_count + j]; @@ -445,7 +448,7 @@ inline std::shared_ptr readDataArr(gdalImage& imgds, long start_row, long sta // delete[] temp; // } // else if (gdal_datatype == GDT_Int64) { - // long* temp = new long[rows_count * cols_count]; + // long* temp = new long[pixel_count64]; // demBand->RasterIO(GF_Read, start_col, start_row, cols_count, rows_count, temp, cols_count, //rows_count, gdal_datatype, 0, 0); for (int i = 0; i < rows_count; i++) { for (int j = 0; j < //cols_count; j++) { datamatrix(i, j) = temp[i * cols_count + j]; @@ -454,15 +457,15 @@ inline std::shared_ptr readDataArr(gdalImage& imgds, long start_row, long sta // delete[] temp; // } else if (gdal_datatype == GDT_Float32) { - float* temp = new float[rows_count * cols_count]; + float* temp = new float[pixel_count64]; demBand->RasterIO(GF_Read, start_col, start_row, cols_count, rows_count, temp, cols_count, rows_count, gdal_datatype, 0, 0); - result = std::shared_ptr(new T[rows_count * cols_count], delArrPtr); + result = std::shared_ptr(new T[pixel_count64], delArrPtr); if (method == GDALREADARRCOPYMETHOD::MEMCPYMETHOD) { - std::memcpy(result.get(), temp, rows_count * cols_count); + std::memcpy(result.get(), temp, pixel_count64); } else if (method == GDALREADARRCOPYMETHOD::VARIABLEMETHOD) { - long count = rows_count * cols_count; + long count = pixel_count64; for (long i = 0; i < count; i++) { result.get()[i] = T(temp[i]); } @@ -470,15 +473,15 @@ inline std::shared_ptr readDataArr(gdalImage& imgds, long start_row, long sta delete[] temp; } else if (gdal_datatype == GDT_Float64) { - double* temp = new double[rows_count * cols_count]; + double* temp = new double[pixel_count64]; demBand->RasterIO(GF_Read, start_col, start_row, cols_count, rows_count, temp, cols_count, rows_count, gdal_datatype, 0, 0); - result = std::shared_ptr(new T[rows_count * cols_count], delArrPtr); + result = std::shared_ptr(new T[pixel_count64], delArrPtr); if (method == GDALREADARRCOPYMETHOD::MEMCPYMETHOD) { - std::memcpy(result.get(), temp, rows_count * cols_count); + std::memcpy(result.get(), temp, pixel_count64); } else if (method == GDALREADARRCOPYMETHOD::VARIABLEMETHOD) { - long count = rows_count * cols_count; + long count = pixel_count64; for (long i = 0; i < count; i++) { result.get()[i] = T(temp[i]); } @@ -487,15 +490,15 @@ inline std::shared_ptr readDataArr(gdalImage& imgds, long start_row, long sta } //else if (gdal_datatype == GDT_CFloat32) { // if (std::is_same>::value || std::is_same>::value) { - // float* temp = new float[rows_count * cols_count * 2]; + // float* temp = new float[pixel_count64 * 2]; // demBand->RasterIO(GF_Read, start_col, start_row, cols_count, rows_count, temp, cols_count, rows_count, gdal_datatype, 0, 0); - // result = std::shared_ptr(new T[rows_count * cols_count], delArrPtr); + // result = std::shared_ptr(new T[pixel_count64], delArrPtr); // if (method == GDALREADARRCOPYMETHOD::MEMCPYMETHOD) { - // std::memcpy(result.get(), temp, rows_count * cols_count); + // std::memcpy(result.get(), temp, pixel_count64); // } // else if (method == GDALREADARRCOPYMETHOD::VARIABLEMETHOD) { - // long count = rows_count * cols_count; + // long count = pixel_count64; // for (long i = 0; i < count; i++) { // result.get()[i] = T(temp[i * 2], // temp[i * 2 + 1]); @@ -509,15 +512,15 @@ inline std::shared_ptr readDataArr(gdalImage& imgds, long start_row, long sta //} //else if (gdal_datatype == GDT_CFloat64 ) { // if (std::is_same>::value || std::is_same>::value) { - // double* temp = new double[rows_count * cols_count * 2]; + // double* temp = new double[pixel_count64 * 2]; // demBand->RasterIO(GF_Read, start_col, start_row, cols_count, rows_count, temp, cols_count, rows_count, gdal_datatype, 0, 0); - // result = std::shared_ptr(new T[rows_count * cols_count], delArrPtr); + // result = std::shared_ptr(new T[pixel_count64], delArrPtr); // if (method == GDALREADARRCOPYMETHOD::MEMCPYMETHOD) { - // std::memcpy(result.get(), temp, rows_count * cols_count); + // std::memcpy(result.get(), temp, pixel_count64); // } // else if (method == GDALREADARRCOPYMETHOD::VARIABLEMETHOD) { - // long count = rows_count * cols_count; + // long count = pixel_count64; // for (long i = 0; i < count; i++) { // result.get()[i] = T(temp[i * 2], // temp[i * 2 + 1]); @@ -556,19 +559,21 @@ inline std::shared_ptr readDataArrComplex(gdalImageComplex& imgds, long start rows_count = start_row + rows_count <= imgds.height ? rows_count : imgds.height - start_row; cols_count = start_col + cols_count <= imgds.width ? cols_count : imgds.width - start_col; + int64_t pixel_count64 = static_cast(rows_count) * static_cast(cols_count); + //Eigen::MatrixXd datamatrix(rows_count, cols_count); if (gdal_datatype == GDT_CFloat32) { if (std::is_same>::value || std::is_same>::value) { - float* temp = new float[rows_count * cols_count * 2]; + float* temp = new float[pixel_count64 * 2]; demBand->RasterIO(GF_Read, start_col, start_row, cols_count, rows_count, temp, cols_count, rows_count, gdal_datatype, 0, 0); - result = std::shared_ptr(new T[rows_count * cols_count], delArrPtr); + result = std::shared_ptr(new T[pixel_count64], delArrPtr); if (method == GDALREADARRCOPYMETHOD::MEMCPYMETHOD) { - std::memcpy(result.get(), temp, rows_count * cols_count); + std::memcpy(result.get(), temp, pixel_count64); } else if (method == GDALREADARRCOPYMETHOD::VARIABLEMETHOD) { - long count = rows_count * cols_count; + long count = pixel_count64; for (long i = 0; i < count; i++) { result.get()[i] = T(temp[i * 2], temp[i * 2 + 1]); @@ -582,15 +587,15 @@ inline std::shared_ptr readDataArrComplex(gdalImageComplex& imgds, long start } else if (gdal_datatype == GDT_CFloat64) { if (std::is_same>::value || std::is_same>::value) { - double* temp = new double[rows_count * cols_count * 2]; + double* temp = new double[pixel_count64 * 2]; demBand->RasterIO(GF_Read, start_col, start_row, cols_count, rows_count, temp, cols_count, rows_count, gdal_datatype, 0, 0); - result = std::shared_ptr(new T[rows_count * cols_count], delArrPtr); + result = std::shared_ptr(new T[pixel_count64], delArrPtr); if (method == GDALREADARRCOPYMETHOD::MEMCPYMETHOD) { - std::memcpy(result.get(), temp, rows_count * cols_count); + std::memcpy(result.get(), temp, pixel_count64); } else if (method == GDALREADARRCOPYMETHOD::VARIABLEMETHOD) { - long count = rows_count * cols_count; + long count = pixel_count64; for (long i = 0; i < count; i++) { result.get()[i] = T(temp[i * 2], temp[i * 2 + 1]); diff --git a/GF3CalibrationAndOrthLib/GF3CalibrationAndOrthLib.cpp b/GF3CalibrationAndOrthLib/GF3CalibrationAndOrthLib.cpp index 1905178..3377076 100644 --- a/GF3CalibrationAndOrthLib/GF3CalibrationAndOrthLib.cpp +++ b/GF3CalibrationAndOrthLib/GF3CalibrationAndOrthLib.cpp @@ -1,2 +1,123 @@ #include "GF3CalibrationAndOrthLib.h" +#include "BaseConstVariable.h" +#include "RasterToolBase.h" +#include "ImageOperatorBase.h" +#include +#include +#include +#include "FileOperator.h" +#include "GPUBaseTool.h" +#include "GPUTool.cuh" + +GF3CALIBRATIONANDORTHLIB_EXPORT int GF3_Sigma0_HH2VV(QString in_SigmaHHRasterPath, QString in_IncidencAngleRasterPath, QString out_SigmaVVRasterPath) +{ + // step 1 检查输入合法性 + { + // 1. 检查输入合法性 + if (in_SigmaHHRasterPath.isEmpty() || in_IncidencAngleRasterPath.isEmpty() || out_SigmaVVRasterPath.isEmpty()) + { + qDebug() << "Input or output file path is empty."; + throw std::invalid_argument("Input or output file path is empty."); + return 1; // 代表输入或输出文件路径为空 + } + // 2. 检查输入影像是否存在 + if (!QFile::exists(in_SigmaHHRasterPath) || !QFile::exists(in_IncidencAngleRasterPath)) + { + qDebug() << "Input file does not exist."; + throw std::runtime_error("Input file does not exist."); + return 2; // 代表输入文件不存在 + } + gdalImage inSigmaHHRaster(in_SigmaHHRasterPath); + gdalImage inIncidencAngleRaster(in_IncidencAngleRasterPath); + // 3. 检查输入影像大小是否一致 + if (inSigmaHHRaster.width != inIncidencAngleRaster.width || inSigmaHHRaster.height != inIncidencAngleRaster.height) + { + qDebug() << "Input images have different sizes."; + throw std::runtime_error("Input images have different sizes."); + return 3; // 代表输入影像大小不一致 + } + + // 复制输入文件到输出文件 + copyFile(in_SigmaHHRasterPath, out_SigmaVVRasterPath); + // 4. 检查输出文件是否存在 + if (!QFile::exists(out_SigmaVVRasterPath)) + { + qDebug() << "Output file does not exist."; + throw std::runtime_error("Output file does not exist."); + return 4; // 代表输出文件不存在 + } + } + + // step 2 执行操作 + { + gdalImage inSigmaHHRaster(in_SigmaHHRasterPath); + gdalImage inIncidencAngleRaster(in_IncidencAngleRasterPath); + gdalImage outSigmaVVRaster(out_SigmaVVRasterPath); + + long width = inSigmaHHRaster.width; + long height = inSigmaHHRaster.height; + + int64_t pixel_count64 = static_cast(height) * static_cast(width); + // 读取输入影像数据 + std::shared_ptr sigmaHHRasterData = readDataArr(inSigmaHHRaster, 0, 0, height, width, 1, GDALREADARRCOPYMETHOD::VARIABLEMETHOD); + std::shared_ptr incidenceAngleData = readDataArr(inIncidencAngleRaster, 0, 0, height, width, 1, GDALREADARRCOPYMETHOD::VARIABLEMETHOD); + std::shared_ptr outSigmaVVRasterData = readDataArr(outSigmaVVRaster, 0, 0, height, width, 1, GDALREADARRCOPYMETHOD::VARIABLEMETHOD); + + // 分块执行数据 + int64_t block_count = Memory1MB / 8 * 500; // 500M内存 作为数据分块大小 + + // 创建交换内存 + + //double* sigmaHHRasterDataPtr_H = (double*)mallocCUDAHost(sizeof(double) * block_count); // 主机内存 + double* sigmaHHRasterDataPtr_D = (double*)mallocCUDADevice(sizeof(double) * block_count);// 设备内存 + + //double* incidenceAngleDataPtr_H = (double*)mallocCUDAHost(sizeof(double) * block_count); // 主机内存 + double* incidenceAngleDataPtr_D = (double*)mallocCUDADevice(sizeof(double) * block_count);// 设备内存 + + //double* outSigmaVVRasterDataPtr_H = (double*)mallocCUDAHost(sizeof(double) * block_count); // 主机内存 + double* outSigmaVVRasterDataPtr_D = (double*)mallocCUDADevice(sizeof(double) * block_count);// 设备内存 + + for (int64_t i = 0; i < pixel_count64; i += block_count) + { + int64_t block_size = (i + block_count) > pixel_count64 ? pixel_count64 - i : block_count; + //采用内存拷贝 从主内存到分块内存 + HostToDevice(sigmaHHRasterDataPtr_D, sigmaHHRasterData.get() + i, sizeof(double) * block_size); + HostToDevice(incidenceAngleDataPtr_D, incidenceAngleData.get() + i, sizeof(double) * block_size); + // 执行GPU计算 + + + + + + // 从设备内存到主机内存 + DeviceToHost(outSigmaVVRasterDataPtr_D, outSigmaVVRasterData.get() + i, sizeof(double) * block_size); + + } + + + + + // 强制释放内存 + + FreeCUDADevice(sigmaHHRasterDataPtr_D); + //FreeCUDAHost(sigmaHHRasterDataPtr_H); + FreeCUDADevice(incidenceAngleDataPtr_D); + //FreeCUDAHost(incidenceAngleDataPtr_H); + FreeCUDADevice(outSigmaVVRasterDataPtr_D); + //FreeCUDAHost(outSigmaVVRasterDataPtr_H); + + sigmaHHRasterData.reset(); + incidenceAngleData.reset(); + outSigmaVVRasterData.reset(); + + } + + + + + + + + return -1;// 代表执行成功 +} diff --git a/GF3CalibrationAndOrthLib/GF3CalibrationAndOrthLib.h b/GF3CalibrationAndOrthLib/GF3CalibrationAndOrthLib.h index a755500..578f899 100644 --- a/GF3CalibrationAndOrthLib/GF3CalibrationAndOrthLib.h +++ b/GF3CalibrationAndOrthLib/GF3CalibrationAndOrthLib.h @@ -15,7 +15,14 @@ - +/// +/// 将HH极化的雷达散射系数转换为VV极化的雷达散射系数 +/// +/// +/// +/// +/// +/// extern "C" GF3CALIBRATIONANDORTHLIB_EXPORT int GF3_Sigma0_HH2VV(QString in_SigmaHHRasterPath,QString in_IncidencAngleRasterPath,QString out_SigmaVVRasterPath); diff --git a/GF3CalibrationAndOrthLib/GF3CalibrationGeoCodingFunCUDA.cu b/GF3CalibrationAndOrthLib/GF3CalibrationGeoCodingFunCUDA.cu index 671e20d..06c8b1e 100644 --- a/GF3CalibrationAndOrthLib/GF3CalibrationGeoCodingFunCUDA.cu +++ b/GF3CalibrationAndOrthLib/GF3CalibrationGeoCodingFunCUDA.cu @@ -13,7 +13,7 @@ #include "GPUTool.cuh" -__device__ __host__ float Computrer_polartionConver_rpolf(float inangle, float alpha) +__device__ __host__ float Computrer_polartionConver_rpol_f(float inangle, float alpha) { float tang = 0.0; if (inangle <= 0.0) { @@ -31,7 +31,7 @@ __device__ __host__ float Computrer_polartionConver_rpolf(float inangle, float a return rpol; } -__device__ __host__ double Computrer_polartionConver_rpold(double inangle, double alpha) +__device__ __host__ double Computrer_polartionConver_rpol_d(double inangle, double alpha) { double tang = 0.0; if (inangle <= 0.0) { @@ -48,3 +48,46 @@ __device__ __host__ double Computrer_polartionConver_rpold(double inangle, doubl double rpol = pow((1.0 + 2.0 * tan_val * tan_val), 2) / pow((1.0 + alpha * tan_val * tan_val), 2); return rpol; } + +__host__ __device__ float polartionConver_f(float insig, float inangle, float alpha = 1.0, bool isvv = true) +{ + float rpol = Computrer_polartionConver_rpol_f(inangle, alpha); + + // dB转线性 + float insig_linear = powf(10.0, insig / 10.0); + float osig = 0.0; + if (isvv) { // C: VV->HH + osig = insig_linear * rpol; + // 返回线性值 + return osig; + } + else { // L: HH->VV + osig = insig_linear / rpol; + // 返回dB值 + return 10.0 * log10f(osig); + } +} + + + +__host__ __device__ double polartionConver_d(double insig, double inangle, double alpha = 1.0, bool isvv = true) +{ + double rpol = Computrer_polartionConver_rpol_d(inangle, alpha); + + // dB转线性 + double insig_linear = pow(10.0, insig / 10.0); + double osig = 0.0; + if (isvv) { // C: VV->HH + osig = insig_linear * rpol; + // 返回线性值 + return osig; + } else { // L: HH->VV + osig = insig_linear / rpol; + // 返回dB值 + return 10.0 * log10(osig); + } +} + + + + diff --git a/GF3CalibrationAndOrthLib/GF3CalibrationGeoCodingFunCUDA.cuh b/GF3CalibrationAndOrthLib/GF3CalibrationGeoCodingFunCUDA.cuh index f211274..b29228d 100644 --- a/GF3CalibrationAndOrthLib/GF3CalibrationGeoCodingFunCUDA.cuh +++ b/GF3CalibrationAndOrthLib/GF3CalibrationGeoCodingFunCUDA.cuh @@ -16,7 +16,7 @@ /// 入射角 /// 转换参数 /// -extern __device__ __host__ float Computrer_polartionConver_rpolf(float inangle, float alpha = 1); +extern __device__ __host__ float Computrer_polartionConver_rpol_f(float inangle, float alpha = 1); /// /// 计算极化转换系数 (double) @@ -24,7 +24,28 @@ extern __device__ __host__ float Computrer_polartionConver_rpolf(float inangle, /// 入射角 /// 转换参数 /// -extern __device__ __host__ double Computrer_polartionConver_rpold(double inangle, double alpha = 1); +extern __device__ __host__ double Computrer_polartionConver_rpol_d(double inangle, double alpha = 1); + +/// +/// 后向散射系数变换函数 +/// +/// sigma0 (dB) +/// 入射角 (度) +/// 默认1，基尔霍夫模型，0.6为汤姆森 +/// true表示输入为VV，false为HH +/// +extern __host__ __device__ double polartionConver_d(double insig, double inangle, double alpha = 1.0, bool isvv = true); + +/// +/// 后向散射系数变换函数 +/// +/// sigma0 (dB) +/// 入射角 (度) +/// 默认1，基尔霍夫模型，0.6为汤姆森 +/// true表示输入为VV，false为HH +/// +extern __host__ __device__ float polartionConver_f(float insig, float inangle, float alpha = 1.0, bool isvv = true); +