From 9708114cf1a56362686b51e1a700b56c3930a0fc Mon Sep 17 00:00:00 2001 From: chenzh <3045316072@qq.com> Date: Wed, 26 Feb 2025 09:45:43 +0800 Subject: [PATCH] =?UTF-8?q?=E6=8F=90=E5=89=8D=E5=A4=84=E7=90=86IFFT?= MIME-Version: 1.0 Content-Type: text/plain; charset=UTF-8 Content-Transfer-Encoding: 8bit --- ALLRelease/ALLRelease.vcxproj | 1 + .../BaseTool/ImageOperatorBase.cpp | 118 +++++++++++++++++- .../BaseTool/ImageOperatorBase.h | 3 + .../SimulationSAR/TBPImageAlgCls.cpp | 74 +++++++++++ .../SimulationSAR/TBPImageAlgCls.h | 7 ++ 5 files changed, 201 insertions(+), 2 deletions(-) diff --git a/ALLRelease/ALLRelease.vcxproj b/ALLRelease/ALLRelease.vcxproj index 7f1cd19..0552d3d 100644 --- a/ALLRelease/ALLRelease.vcxproj +++ b/ALLRelease/ALLRelease.vcxproj @@ -159,6 +159,7 @@ true stdcpp14 stdc11 + true Console diff --git a/BaseCommonLibrary/BaseTool/ImageOperatorBase.cpp b/BaseCommonLibrary/BaseTool/ImageOperatorBase.cpp index 1cae5e5..f343442 100644 --- a/BaseCommonLibrary/BaseTool/ImageOperatorBase.cpp +++ b/BaseCommonLibrary/BaseTool/ImageOperatorBase.cpp @@ -1822,6 +1822,28 @@ gdalImageComplex CreategdalImageComplex(const QString& img_path, int height, int return result_img; } +gdalImageComplex BASECONSTVARIABLEAPI CreategdalImageComplexNoProj(const QString& img_path, int height, int width, int band_num, bool overwrite) +{ + + // 创建新文件 + omp_lock_t lock; + omp_init_lock(&lock); + omp_set_lock(&lock); + GDALAllRegister(); + CPLSetConfigOption("GDAL_FILENAME_IS_UTF8", "YES"); + + GDALDriver* poDriver = GetGDALDriverManager()->GetDriverByName("ENVI"); + GDALDataset* poDstDS = (poDriver->Create(img_path.toUtf8().constData(), width, height, band_num, GDT_CFloat64, NULL)); + GDALFlushCache((GDALDatasetH)poDstDS); + GDALClose((GDALDatasetH)poDstDS); + //poDstDS.reset(); + omp_unset_lock(&lock); // + omp_destroy_lock(&lock); // + + gdalImageComplex result_img(img_path); + return result_img; +} + gdalImageComplex CreateEchoComplex(const QString& img_path, int height, int width, int band_num) { // 创建图像 @@ -2598,6 +2620,60 @@ void gdalImageComplex::saveImage(Eigen::MatrixXcd data, int start_row, int start omp_destroy_lock(&lock); // } +void gdalImageComplex::saveImage(std::shared_ptr> data, long start_row, long start_col, long rowCount, long colCount, int band_ids) +{ + omp_lock_t lock; + omp_init_lock(&lock); + omp_set_lock(&lock); + if (start_row + rowCount > this->height || start_col + colCount > this->width) { + QString tip = u8"file path: " + this->img_path; + qDebug() << tip; + throw std::exception(tip.toUtf8().constData()); + return; + } + GDALAllRegister(); + CPLSetConfigOption("GDAL_FILENAME_IS_UTF8", "YES"); + GDALDriver* poDriver = GetGDALDriverManager()->GetDriverByName("ENVI"); + GDALDataset* poDstDS = nullptr; + if (exists_test(this->img_path)) { + poDstDS = (GDALDataset*)(GDALOpen(this->img_path.toUtf8().constData(), GA_Update)); + } + else { + poDstDS = poDriver->Create(this->img_path.toUtf8().constData(), this->width, this->height, + this->band_num, GDT_CFloat64, NULL); // 斤拷锟斤拷 + poDstDS->SetProjection(this->projection.toUtf8().constData()); + + double gt_ptr[6]; + for (int i = 0; i < this->gt.rows(); i++) { + for (int j = 0; j < this->gt.cols(); j++) { + gt_ptr[i * 3 + j] = this->gt(i, j); + } + } + poDstDS->SetGeoTransform(gt_ptr); + //delete[] gt_ptr; + } + + double* databuffer = new double[rowCount * colCount * 2]; + for (long i = 0; i < rowCount; i++) { + for (long j = 0; j < colCount; j++) { + databuffer[i * colCount * 2 + j * 2] = data.get()[i*colCount+j].real(); + databuffer[i * colCount * 2 + j * 2 + 1] = data.get()[i * colCount + j].imag(); + } + } + + // poDstDS->RasterIO(GF_Write,start_col, start_row, datacols, datarows, databuffer, datacols, + // datarows, GDT_Float32,band_ids, num,0,0,0); + poDstDS->GetRasterBand(band_ids)->RasterIO(GF_Write, start_col, start_row, colCount, rowCount, + databuffer, colCount, rowCount, GDT_CFloat64, 0, 0); + GDALFlushCache(poDstDS); + delete databuffer; + GDALClose((GDALDatasetH)poDstDS); + // GDALDestroy(); // or, DllMain at DLL_PROCESS_DETACH + omp_unset_lock(&lock); // + omp_destroy_lock(&lock); // + +} + Eigen::MatrixXcd gdalImageComplex::getDataComplex(int start_row, int start_col, int rows_count, int cols_count, int band_ids) { @@ -2617,8 +2693,8 @@ Eigen::MatrixXcd gdalImageComplex::getDataComplex(int start_row, int start_col, poBand = poDataset->GetRasterBand(1); // 读取波段信息，假设是复数类型 - int nXSize = poBand->GetXSize(); - int nYSize = poBand->GetYSize(); + int nXSize = cols_count; poBand->GetXSize(); + int nYSize = rows_count; poBand->GetYSize(); double* databuffer = new double[nXSize * nYSize * 2]; poBand->RasterIO(GF_Read, start_col, start_row, cols_count, rows_count, databuffer, cols_count, @@ -2635,6 +2711,44 @@ Eigen::MatrixXcd gdalImageComplex::getDataComplex(int start_row, int start_col, delete[] databuffer; return rasterData; } + +std::shared_ptr> gdalImageComplex::getDataComplexSharePtr(int start_row, int start_col, int rows_count, int cols_count, int band_ids) +{ + GDALDataset* poDataset; + GDALAllRegister(); + CPLSetConfigOption("GDAL_FILENAME_IS_UTF8", "YES"); + + // 打开TIFF文件 + poDataset = (GDALDataset*)GDALOpen(this->img_path.toUtf8().constData(), GA_ReadOnly); + if(poDataset == nullptr) { + QMessageBox::warning(nullptr, u8"错误", u8"无法打开文件：" + this->img_path); + qDebug() << u8"无法打开文件：" + this->img_path; + } + + // 获取数据集的第一个波段 + GDALRasterBand* poBand; + poBand = poDataset->GetRasterBand(1); + + // 读取波段信息，假设是复数类型 + + + double* databuffer = new double[rows_count * cols_count * 2]; + poBand->RasterIO(GF_Read, start_col, start_row, cols_count, rows_count, databuffer, cols_count, + rows_count, GDT_CFloat64, 0, 0); + GDALClose((GDALDatasetH)poDataset); + + std::shared_ptr> rasterData(new std::complex[rows_count * cols_count], delArrPtr); + + for(size_t i = 0; i < rows_count; i++) { + for(size_t j = 0; j < cols_count; j++) { + rasterData.get()[i*cols_count+ j] = std::complex(databuffer[i * cols_count * 2 + j * 2], + databuffer[i * cols_count * 2 + j * 2 + 1]); + } + } + + delete[] databuffer; + return rasterData; +} void gdalImageComplex::saveImage() { diff --git a/BaseCommonLibrary/BaseTool/ImageOperatorBase.h b/BaseCommonLibrary/BaseTool/ImageOperatorBase.h index 1b5a223..4c63230 100644 --- a/BaseCommonLibrary/BaseTool/ImageOperatorBase.h +++ b/BaseCommonLibrary/BaseTool/ImageOperatorBase.h @@ -228,7 +228,9 @@ public: // 方法 ~gdalImageComplex(); void setData(Eigen::MatrixXcd); void saveImage(Eigen::MatrixXcd data, int start_row, int start_col, int band_ids); + void saveImage(std::shared_ptr> data, long start_row, long start_col, long rowCount, long colCount, int band_ids); Eigen::MatrixXcd getDataComplex(int start_row, int start_col, int rows_count, int cols_count, int band_ids); + std::shared_ptr> getDataComplexSharePtr(int start_row, int start_col, int rows_count, int cols_count, int band_ids); void saveImage() override; void savePreViewImage(); @@ -243,6 +245,7 @@ gdalImage BASECONSTVARIABLEAPI CreategdalImage(const QString& img_path, int gdalImage BASECONSTVARIABLEAPI CreategdalImage(const QString& img_path, int height, int width, int band_num, Eigen::MatrixXd gt, long espgcode, GDALDataType eType = GDT_Float32, bool need_gt = true, bool overwrite = false, bool isENVI = false); gdalImageComplex BASECONSTVARIABLEAPI CreategdalImageComplex(const QString& img_path, int height, int width, int band_num, Eigen::MatrixXd gt, QString projection, bool need_gt = true, bool overwrite = false); +gdalImageComplex BASECONSTVARIABLEAPI CreategdalImageComplexNoProj(const QString& img_path, int height, int width, int band_num, bool overwrite = true); gdalImageComplex BASECONSTVARIABLEAPI CreateEchoComplex(const QString& img_path, int height, int width, int band_num); diff --git a/Toolbox/SimulationSARTool/SimulationSAR/TBPImageAlgCls.cpp b/Toolbox/SimulationSARTool/SimulationSAR/TBPImageAlgCls.cpp index 3041fbe..f7f5078 100644 --- a/Toolbox/SimulationSARTool/SimulationSAR/TBPImageAlgCls.cpp +++ b/Toolbox/SimulationSARTool/SimulationSAR/TBPImageAlgCls.cpp @@ -203,6 +203,12 @@ ErrorCode TBPImageAlgCls::Process(long num_thread) this->L1ds->saveImageRaster(imageRaster, startline,templine); } qDebug() << u8"ʼ"; + qDebug() << u8"Ƶز-> ʱز"; + this->TimeEchoDataPath = JoinPath(this->L1ds->getoutFolderPath(), this->L0ds->getSimulationTaskName() + "_Timeecho.bin"); + this->EchoFreqToTime(); + + qDebug() << u8"Ƶز-> ʱز "; + if (GPURUN) { return this->ProcessGPU(); } @@ -361,6 +367,74 @@ ErrorCode TBPImageAlgCls::ProcessGPU() return ErrorCode::SUCCESS; } +void TBPImageAlgCls::EchoFreqToTime( ) +{ + // ȡ + + long PRFCount = this->L0ds->getPluseCount(); + long inColCount = this->L0ds->getPlusePoints(); + long outColCount = nextpow2(inColCount); + this->TimeEchoRowCount = PRFCount; + this->TimeEchoColCount = outColCount; + qDebug() << "IFFT : " << this->TimeEchoDataPath; + qDebug() << "PRF Count:\t" << PRFCount; + qDebug() << "inColCount:\t" << inColCount; + qDebug() << "outColCount:\t" << outColCount; + // ļ + gdalImageComplex outTimeEchoImg = CreategdalImageComplexNoProj(this->TimeEchoDataPath,this->TimeEchoRowCount,this->TimeEchoColCount,1); + + // ֿ + long echoBlockline = Memory1GB / 8 / 2 / outColCount * 1; //1GB + echoBlockline = echoBlockline < 1 ? 1 : echoBlockline; + + + long startechoid = 0; + for (long startechoid = 0; startechoid < PRFCount; startechoid = startechoid + echoBlockline) { + + long tempechoBlockline = echoBlockline; + if (startechoid + tempechoBlockline >= PRFCount) { + tempechoBlockline = PRFCount - startechoid; + } + std::shared_ptr> echoArr = this->L0ds->getEchoArr(startechoid, tempechoBlockline); + std::shared_ptr> IFFTArr = outTimeEchoImg.getDataComplexSharePtr(startechoid, 0, tempechoBlockline, outColCount, 1); + + std::shared_ptr host_echoArr((cuComplex*)mallocCUDAHost(sizeof(cuComplex)* tempechoBlockline * inColCount), FreeCUDAHost); + std::shared_ptr host_IFFTechoArr((cuComplex*)mallocCUDAHost(sizeof(cuComplex)* tempechoBlockline * outColCount), FreeCUDAHost); + +#pragma omp parallel for + for (long ii = 0; ii < tempechoBlockline * inColCount; ii++) { + host_echoArr.get()[ii] = make_cuComplex(echoArr.get()[ii].real(), echoArr.get()[ii].imag()); + } +#pragma omp parallel for + for (long ii = 0; ii < tempechoBlockline * outColCount; ii++) { + host_IFFTechoArr.get()[ii] = make_cuComplex(0,0); + } + + std::shared_ptr device_echoArr((cuComplex*)mallocCUDADevice(sizeof(cuComplex) * tempechoBlockline * inColCount), FreeCUDADevice); + std::shared_ptr device_IFFTechoArr((cuComplex*)mallocCUDADevice(sizeof(cuComplex) * tempechoBlockline * outColCount), FreeCUDADevice); + + HostToDevice(host_echoArr.get(), device_echoArr.get(), sizeof(cuComplex) * tempechoBlockline * inColCount); + HostToDevice(host_IFFTechoArr.get(), device_IFFTechoArr.get(), sizeof(cuComplex) * tempechoBlockline * outColCount); + CUDAIFFT(device_echoArr.get(), device_IFFTechoArr.get(), tempechoBlockline, inColCount, outColCount); + + DeviceToHost(host_IFFTechoArr.get(), device_IFFTechoArr.get(), sizeof(cuComplex) * tempechoBlockline * outColCount); + +#pragma omp parallel for + for (long ii = 0; ii < tempechoBlockline * outColCount; ii++) { + IFFTArr.get()[ii] = std::complex(host_IFFTechoArr.get()[ii].x, host_IFFTechoArr.get()[ii].y); + } + + outTimeEchoImg.saveImage(IFFTArr, startechoid, 0, tempechoBlockline, outColCount, 1); + + qDebug() << QString(" image block PRF:[%1] \t").arg((startechoid + tempechoBlockline) * 100.0 / PRFCount) + << startechoid << "\t-\t" << startechoid + tempechoBlockline; + } + + return; + + +} + diff --git a/Toolbox/SimulationSARTool/SimulationSAR/TBPImageAlgCls.h b/Toolbox/SimulationSARTool/SimulationSAR/TBPImageAlgCls.h index 9f8f753..a0f9332 100644 --- a/Toolbox/SimulationSARTool/SimulationSAR/TBPImageAlgCls.h +++ b/Toolbox/SimulationSARTool/SimulationSAR/TBPImageAlgCls.h @@ -57,6 +57,13 @@ private: //ErrorCode ProcessCPU(long num_thread); ErrorCode ProcessGPU(); +private://ʱԱ + QString TimeEchoDataPath; + long TimeEchoRowCount; + long TimeEchoColCount; + + void EchoFreqToTime( ); + };