From d45b3975608e479fb07cc387fc8d4784779777b2 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?=E9=99=88=E5=A2=9E=E8=BE=89?= <3045316072@qq.com> Date: Fri, 18 Apr 2025 16:38:06 +0800 Subject: [PATCH] =?UTF-8?q?=E4=BF=AE=E6=AD=A3=E5=A4=9A=E6=99=AE=E5=8B=92?= =?UTF-8?q?=E8=AE=A1=E7=AE=97=E5=85=AC=E5=BC=8F=EF=BC=88=E5=8E=BB=E6=8E=89?= =?UTF-8?q?=E5=9C=B0=E9=80=9F=EF=BC=89?= MIME-Version: 1.0 Content-Type: text/plain; charset=UTF-8 Content-Transfer-Encoding: 8bit --- .../BaseTool/SARSimulationImageL1.cpp | 17 +- .../BaseTool/SARSimulationImageL1.h | 2 + BaseCommonLibrary/ImageOperatorFuntion.cpp | 89 +++++------ Toolbox/BaseToolbox/BaseToolbox.vcxproj | 1 + .../GF3CalibrationAndGeocodingClass.cpp | 147 ++++++++++++------ .../BaseToolbox/QOrthSlrRaster.cpp | 8 +- 6 files changed, 163 insertions(+), 101 deletions(-) diff --git a/BaseCommonLibrary/BaseTool/SARSimulationImageL1.cpp b/BaseCommonLibrary/BaseTool/SARSimulationImageL1.cpp index 8bc4245..2fbfbb8 100644 --- a/BaseCommonLibrary/BaseTool/SARSimulationImageL1.cpp +++ b/BaseCommonLibrary/BaseTool/SARSimulationImageL1.cpp @@ -251,6 +251,7 @@ void SARSimulationImageL1Dataset::saveToXml() xmlWriter.writeTextElement("Rref", QString::number(this->Rref, 'e', 18)); xmlWriter.writeTextElement("CenterFreq", QString::number(this->centerFreq, 'e', 18)); xmlWriter.writeTextElement("Fs", QString::number(this->Fs, 'e', 18)); + xmlWriter.writeTextElement("PRF", QString::number(this->prf, 'e', 18)); xmlWriter.writeTextElement("CenterAngle", QString::number(this->CenterAngle, 'e', 18)); xmlWriter.writeTextElement("LookSide", this->LookSide); @@ -371,6 +372,9 @@ ErrorCode SARSimulationImageL1Dataset::loadFromXml() else if (xmlReader.name() == "Fs") { this->Fs = xmlReader.readElementText().toDouble(); } + else if (xmlReader.name() == "PRF") { + this->prf = xmlReader.readElementText().toDouble(); + } else if(xmlReader.name() == "ImageStartTime"){ this->startImageTime = xmlReader.readElementText().toDouble(); } @@ -873,6 +877,7 @@ QVector SARSimulationImageL1Dataset::getDopplerParams() result[2] = d2; result[3] = d3; result[4] = d4; + return result; } @@ -888,11 +893,13 @@ void SARSimulationImageL1Dataset::setDopplerParams(double id0, double id1, doubl QVector SARSimulationImageL1Dataset::getDopplerCenterCoff() { QVector result(5); - result[0]=r0; - result[1]=r1; - result[2]=r2; - result[3]=r3; - result[4]=r4; + result[0] = r0; + result[1] = r1; + result[2] = r2; + result[3] = r3; + result[4] = r4; + + return result; } diff --git a/BaseCommonLibrary/BaseTool/SARSimulationImageL1.h b/BaseCommonLibrary/BaseTool/SARSimulationImageL1.h index 3d85682..6d8325a 100644 --- a/BaseCommonLibrary/BaseTool/SARSimulationImageL1.h +++ b/BaseCommonLibrary/BaseTool/SARSimulationImageL1.h @@ -143,9 +143,11 @@ public: double getDopplerParametersReferenceTime(); void setDopplerParametersReferenceTime(double v); + // 多普勒参数 QVector getDopplerParams(); void setDopplerParams(double d0, double d1, double d2, double d3, double d4); + // 多普勒中心系数 QVector getDopplerCenterCoff(); void setDopplerCenterCoff(double r0, double r1, double r2, double r3, double r4); diff --git a/BaseCommonLibrary/ImageOperatorFuntion.cpp b/BaseCommonLibrary/ImageOperatorFuntion.cpp index ae75a14..a7ec078 100644 --- a/BaseCommonLibrary/ImageOperatorFuntion.cpp +++ b/BaseCommonLibrary/ImageOperatorFuntion.cpp @@ -633,8 +633,10 @@ void cropRasterByLatLon(const char* inputFile, const char* outputFile, double mi qDebug() << "Raster cropped and saved to: " << outputFile; // 娓呯悊 - GDALClose(poDataset); - GDALClose(poOutDataset); + + GDALClose((GDALDatasetH)(GDALDatasetH)poDataset); + GDALClose((GDALDatasetH)(GDALDatasetH)poOutDataset); + GDALDestroy(); // or, DllMain at DLL_PROCESS_DETACH } ErrorCode transformCoordinate(double x, double y, int sourceEPSG, int targetEPSG, Point2& p) { @@ -763,8 +765,8 @@ void resampleRaster(const char* inputRaster, const char* outputRaster, double ta int nYSize = poDataset->GetRasterYSize(); // 璁$畻鐩爣鏍呮牸鐨勫昂瀵 - double targetXSize = (adfGeoTransform[1] * nXSize) / targetPixelSizeX; - double targetYSize = std::abs((adfGeoTransform[5] * nYSize) / targetPixelSizeY); + double targetXSize = ceil((adfGeoTransform[1] * nXSize) / targetPixelSizeX); + double targetYSize = ceil(std::abs((adfGeoTransform[5] * nYSize) / targetPixelSizeY)); // 鍒涘缓鐩爣鏁版嵁闆嗭紙杈撳嚭鏍呮牸锛 GDALDriver* poDriver = GetGDALDriverManager()->GetDriverByName("GTiff"); @@ -774,38 +776,13 @@ void resampleRaster(const char* inputRaster, const char* outputRaster, double ta return; } - // 鍒涘缓杈撳嚭鏁版嵁闆 - GDALDataset* poOutDataset = poDriver->Create(outputRaster, (int)targetXSize, (int)targetYSize, poDataset->GetRasterCount(), GDT_Float32, nullptr); - if (poOutDataset == nullptr) { - qDebug() << "Failed to create output raster." ; - GDALClose(poDataset); - return; - } - // 璁剧疆杈撳嚭鏁版嵁闆嗙殑鍦扮悊鍙樻崲锛堝潗鏍囩郴锛 double targetGeoTransform[6] = { adfGeoTransform[0], targetPixelSizeX, 0, adfGeoTransform[3], 0, -targetPixelSizeY }; - poOutDataset->SetGeoTransform(targetGeoTransform); - - // 璁剧疆杈撳嚭鏁版嵁闆嗙殑鎶曞奖淇℃伅 - poOutDataset->SetProjection(poDataset->GetProjectionRef()); - - // 杩涜閲嶉噰鏍 - for (int i = 0; i < poDataset->GetRasterCount(); i++) { - GDALRasterBand* poBand = poDataset->GetRasterBand(i + 1); - GDALRasterBand* poOutBand = poOutDataset->GetRasterBand(i + 1); - - // 浣跨敤GDAL鐨勯噸閲囨牱鏂规硶锛岄夋嫨涓涓傚綋鐨勯噸閲囨牱鏂瑰紡 - poOutBand->RasterIO(GF_Write, 0, 0, (int)targetXSize, (int)targetYSize, - nullptr, (int)targetXSize, (int)targetYSize, - poBand->GetRasterDataType(), 0, 0, nullptr); - } - - // 鍏抽棴鏁版嵁闆 - GDALClose(poDataset); - GDALClose(poOutDataset); - + GDALClose((GDALDatasetH)(GDALDatasetH)poDataset); + GDALDestroy(); // or, DllMain at DLL_PROCESS_DETACH + ResampleGDAL(inputRaster, outputRaster, targetGeoTransform, targetXSize, targetYSize,GRA_Bilinear); qDebug() << "Resampling completed."; } @@ -1438,23 +1415,30 @@ void clipRaster(QString inRasterPath, QString outRasterPath, long minRow, long m ErrorCode DEM2XYZRasterAndSlopRaster(QString dempath, QString demxyzpath, QString demsloperPath) { - gdalImage demds(dempath); - gdalImage demxyz = CreategdalImage(demxyzpath, demds.height, demds.width, 3, demds.gt, demds.projection, true, true);// X,Y,Z + QString DEMPath = dempath; + QString XYZPath = demxyzpath; + QString SLOPERPath = demsloperPath; + + + + gdalImage demds(DEMPath); + gdalImage demxyz = CreategdalImageDouble(XYZPath, demds.height, demds.width, 3, demds.gt, demds.projection, true, true, false);// X,Y,Z // 鍒嗗潡璁$畻骞惰浆鎹负XYZ - Eigen::MatrixXd demArr = demds.getData(0, 0, demds.height, demds.width, 1); - Eigen::MatrixXd demR = demArr; - Landpoint LandP{ 0,0,0 }; - Point3 GERpoint{ 0,0,0 }; + //Eigen::MatrixXd demArr = demds.getData(0, 0, demds.height, demds.width, 1); + //Eigen::MatrixXd demR = demArr; + double R = 0; double dem_row = 0, dem_col = 0, dem_alt = 0; - long line_invert = 1000; + long line_invert = Memory1MB / 8.0 / demds.width * 1000; + + - double rowidx = 0; - double colidx = 0; for (int max_rows_ids = 0; max_rows_ids < demds.height; max_rows_ids = max_rows_ids + line_invert) { + + Eigen::MatrixXd demdata = demds.getData(max_rows_ids, 0, line_invert, demds.width, 1); Eigen::MatrixXd xyzdata_x = demdata.array() * 0; Eigen::MatrixXd xyzdata_y = demdata.array() * 0; @@ -1463,7 +1447,12 @@ ErrorCode DEM2XYZRasterAndSlopRaster(QString dempath, QString demxyzpath, QStri int datarows = demdata.rows(); int datacols = demdata.cols(); +#pragma omp parallel for for (int i = 0; i < datarows; i++) { + Landpoint LandP{ 0,0,0 }; + Point3 GERpoint{ 0,0,0 }; + double rowidx = 0; + double colidx = 0; for (int j = 0; j < datacols; j++) { rowidx = i + max_rows_ids; colidx = j; @@ -1477,19 +1466,23 @@ ErrorCode DEM2XYZRasterAndSlopRaster(QString dempath, QString demxyzpath, QStri demxyz.saveImage(xyzdata_x, max_rows_ids, 0, 1); demxyz.saveImage(xyzdata_y, max_rows_ids, 0, 2); demxyz.saveImage(xyzdata_z, max_rows_ids, 0, 3); + } // 璁$畻鍧″悜瑙 - gdalImage demsloperxyz = CreategdalImage(demsloperPath, demds.height, demds.width, 4, demds.gt, demds.projection, true, true);// X,Y,Z,cosangle + + + + gdalImage demsloperxyz = CreategdalImageDouble(SLOPERPath, demds.height, demds.width, 4, demds.gt, demds.projection, true, true, false);// X,Y,Z,cosangle line_invert = 1000; long start_ids = 0; long dem_rows = 0, dem_cols = 0; - for (start_ids = 1; start_ids < demds.height; start_ids = start_ids + line_invert) { Eigen::MatrixXd demdata = demds.getData(start_ids - 1, 0, line_invert + 2, demxyz.width, 1); long startlineid = start_ids; + Eigen::MatrixXd demsloper_x = demsloperxyz.getData(start_ids - 1, 0, line_invert + 2, demxyz.width, 1); Eigen::MatrixXd demsloper_y = demsloperxyz.getData(start_ids - 1, 0, line_invert + 2, demxyz.width, 2); Eigen::MatrixXd demsloper_z = demsloperxyz.getData(start_ids - 1, 0, line_invert + 2, demxyz.width, 3); @@ -1498,13 +1491,13 @@ ErrorCode DEM2XYZRasterAndSlopRaster(QString dempath, QString demxyzpath, QStri Landpoint p0, p1, p2, p3, p4, pslopeVector, pp; Vector3D slopeVector; - dem_rows = demsloper_y.rows(); - dem_cols = demsloper_y.cols(); + dem_rows = demsloper_x.rows(); + dem_cols = demsloper_x.cols(); double sloperAngle = 0; Vector3D Zaxis = { 0,0,1 }; double rowidx = 0, colidx = 0; - + //#pragma omp parallel for for (long i = 1; i < dem_rows - 1; i++) { for (long j = 1; j < dem_cols - 1; j++) { rowidx = i + startlineid; @@ -1527,15 +1520,19 @@ ErrorCode DEM2XYZRasterAndSlopRaster(QString dempath, QString demxyzpath, QStri demsloper_y(i, j) = slopeVector.y; demsloper_z(i, j) = slopeVector.z; demsloper_angle(i, j) = sloperAngle; + } } + demsloperxyz.saveImage(demsloper_x, start_ids - 1, 0, 1); demsloperxyz.saveImage(demsloper_y, start_ids - 1, 0, 2); demsloperxyz.saveImage(demsloper_z, start_ids - 1, 0, 3); demsloperxyz.saveImage(demsloper_angle, start_ids - 1, 0, 4); + } + return ErrorCode::SUCCESS; } diff --git a/Toolbox/BaseToolbox/BaseToolbox.vcxproj b/Toolbox/BaseToolbox/BaseToolbox.vcxproj index 699771e..7fe5c9a 100644 --- a/Toolbox/BaseToolbox/BaseToolbox.vcxproj +++ b/Toolbox/BaseToolbox/BaseToolbox.vcxproj @@ -122,6 +122,7 @@ true stdcpp14 stdc11 + Disabled diff --git a/Toolbox/BaseToolbox/BaseToolbox/GF3CalibrationAndGeocodingClass.cpp b/Toolbox/BaseToolbox/BaseToolbox/GF3CalibrationAndGeocodingClass.cpp index dc3c8f9..a78cae3 100644 --- a/Toolbox/BaseToolbox/BaseToolbox/GF3CalibrationAndGeocodingClass.cpp +++ b/Toolbox/BaseToolbox/BaseToolbox/GF3CalibrationAndGeocodingClass.cpp @@ -212,7 +212,7 @@ ErrorCode ImportGF3L1AProcess(QString inMetaxmlPath, QString outWorkDirPath) ErrorCode RD_PSTN(double& refrange, double& lamda, double& timeR, double& R, double& tx, double& ty, double& tz, double& slopex, double& slopey, double& slopez, GF3PolyfitSatelliteOribtModel& polyfitmodel, SatelliteOribtNode& node, double& d0, double& d1, double& d2, double& d3, double& d4) { - double dt = 1e-6; + double dt = 1e-4; double inct = 0; bool antfalg = false; double timeR2 = 0; @@ -231,21 +231,21 @@ ErrorCode RD_PSTN(double& refrange, double& lamda, double& timeR, double& R, dou for (int i = 0; i < 50; i++) { polyfitmodel.getSatelliteOribtNode(timeR, node, antfalg); R1 = std::sqrt(std::pow(node.Px - tx, 2) + std::pow(node.Py - ty, 2) + std::pow(node.Pz - tz, 2)); - dplerTheory1 = (-2 / lamda) * (((node.Px - tx) * (node.Vx + EARTHWE * ty) + (node.Py - ty) * (node.Vy - EARTHWE * tz) + (node.Pz - tz) * (node.Vz - 0)) / R1); - dplerR = (R1 - refrange) * 1e6 / LIGHTSPEED; + dplerTheory1 = (-2 / lamda / R1) * ((node.Px - tx) * (node.Vx /*+ EARTHWE * ty*/) + (node.Py - ty) * (node.Vy /*- EARTHWE * tz*/) + (node.Pz - tz) * (node.Vz - 0)); + dplerR = (R1 - refrange) * 1e6 / LIGHTSPEED; // GF3 计算公式 dplerNumber1 = d0 + dplerR * d1 + std::pow(dplerR, 2) * d2 + std::pow(dplerR, 3) * d3 + std::pow(dplerR, 4) * d4; timeR2 = timeR + dt; polyfitmodel.getSatelliteOribtNode(timeR2, node, antfalg); R2 = std::sqrt(std::pow(node.Px - tx, 2) + std::pow(node.Py - ty, 2) + std::pow(node.Pz - tz, 2)); - dplerTheory2 = (-2 / lamda) * (((node.Px - tx) * (node.Vx + EARTHWE * ty) + (node.Py - ty) * (node.Vy - EARTHWE * tz) + (node.Pz - tz) * (node.Vz - 0)) / R2); - dplerR = (R2 - refrange) * 1e6 / LIGHTSPEED; - dplerNumber2 = d0 + dplerR * d1 + std::pow(dplerR, 2) * d2 + std::pow(dplerR, 3) * d3 + std::pow(dplerR, 4) * d4; - - inct = dt * (dplerTheory2 - dplerNumber1) / (dplerTheory1 - dplerTheory2); - if (std::abs(dplerTheory1 - dplerTheory2) < 1e-9 || std::abs(inct) < dt) { + dplerTheory2 = (-2 / lamda/R2) * ((node.Px - tx) * (node.Vx/* + EARTHWE * ty*/) + (node.Py - ty) * (node.Vy /*- EARTHWE * tz*/) + (node.Pz - tz) * (node.Vz - 0)); + + + inct = (dplerTheory2-dplerNumber1) / (dplerTheory2-dplerTheory1 ); + if (std::abs(dplerNumber1 - dplerTheory2) < 1e-6 || std::abs(inct) < 1.0e-4) { break; } + inct = std::abs(inct) < 10 ?inct*1.0e-2:inct*dt; timeR = timeR - inct; } R = R1; // 斜距 @@ -265,7 +265,7 @@ ErrorCode GF3RDCreateLookTable(QString inxmlPath, QString indemPath, QString out gdalImage rasterRC = CreategdalImage(outlooktablePath, demxyz.height, demxyz.width, 2, demxyz.gt, demxyz.projection, true, true);// X,Y,Z gdalImage localincangleRaster; if (localincAngleFlag) { - localincangleRaster = CreategdalImage(outLocalIncidenceAnglePath, demxyz.height, demxyz.width, 2, demxyz.gt, demxyz.projection, true, true);// X,Y,Z + localincangleRaster = CreategdalImage(outLocalIncidenceAnglePath, demxyz.height, demxyz.width, 1, demxyz.gt, demxyz.projection, true, true);// X,Y,Z } @@ -279,28 +279,48 @@ ErrorCode GF3RDCreateLookTable(QString inxmlPath, QString indemPath, QString out double d3 = dopplers[3]; double d4 = dopplers[4]; - double fs = l1dataset.getFs();// Fs采样 - double prf = (l1dataset.getEndImageTime() - l1dataset.getStartImageTime()) / (l1dataset.getrowCount() - 1);// PRF 采样 + double fs = l1dataset.getFs()*1e6;// Fs采样 + double prf = l1dataset.getPRF();// PRF 采样 double nearRange = l1dataset.getNearRange(); double imagestarttime = l1dataset.getStartImageTime(); + double imageendtime = l1dataset.getEndImageTime(); double refrange = l1dataset.getRefRange(); - double lamda = (LIGHTSPEED*1e-6)/ l1dataset.getCenterFreq(); + double lamda = (LIGHTSPEED*1.0e-9)/ l1dataset.getCenterFreq(); + // 打印参数 + qDebug() << "Fs:\t" << fs; + qDebug() << "prf:\t" << prf; + qDebug() << "nearRange:\t" << nearRange; + qDebug() << "imagestarttime:\t" << imagestarttime; + qDebug() << "imageendtime:\t" << imageendtime; + qDebug() << "refrange:\t" << refrange; + qDebug() << "lamda:\t" << lamda; + //打印多普累参数 + qDebug() << u8"-----------------------------------"; + qDebug() << u8"多普勒参数:\n"; + qDebug() << u8"d0:\t" << d0; + qDebug() << u8"d1:\t" << d1; + qDebug() << u8"d2:\t" << d2; + qDebug() << u8"d3:\t" << d3; + qDebug() << u8"d4:\t" << d4; + qDebug() << u8"-----------------------------------"; // 构建轨道模型 GF3PolyfitSatelliteOribtModel polyfitmodel; QVector < SatelliteAntPos > antposes = l1dataset.getXmlSateAntPos(); - polyfitmodel.setSatelliteOribtStartTime(imagestarttime); + polyfitmodel.setSatelliteOribtStartTime((imagestarttime+imageendtime)/2); for (long i = 0; i < antposes.size(); i++) { - SatelliteOribtNode node; - node.time = antposes[i].time; - node.Px = antposes[i].Px; - node.Py = antposes[i].Py; - node.Pz = antposes[i].Pz; - node.Vx = antposes[i].Vx; - node.Vy = antposes[i].Vy; - node.Vz = antposes[i].Vz; - polyfitmodel.addOribtNode(node); + if (antposes[i].time > imagestarttime - 5 && antposes[i].time < imageendtime + 5) { + SatelliteOribtNode node; + node.time = antposes[i].time; + node.Px = antposes[i].Px; + node.Py = antposes[i].Py; + node.Pz = antposes[i].Pz; + node.Vx = antposes[i].Vx; + node.Vy = antposes[i].Vy; + node.Vz = antposes[i].Vz; + polyfitmodel.addOribtNode(node); + } } polyfitmodel.polyFit(3, false); @@ -312,7 +332,7 @@ ErrorCode GF3RDCreateLookTable(QString inxmlPath, QString indemPath, QString out // 开始计算查找表 //1.计算分块 long cpucore_num = std::thread::hardware_concurrency(); - long blocklineinit = Memory1MB / 8 / cpucore_num / 4 / l1dataset.getcolCount() * 4000; + long blocklineinit = Memory1MB / 8 / cpucore_num / 4 / l1dataset.getcolCount() * 8000; blocklineinit = blocklineinit < 50 ? 50 : blocklineinit; //2.迭代计算 long colcount = rasterRC.width;//l1dataset.getcolCount(); @@ -333,7 +353,7 @@ ErrorCode GF3RDCreateLookTable(QString inxmlPath, QString indemPath, QString out -#pragma omp parallel for num_threads(cpucore_num-1) + for (startRId = 0; startRId < rowcount; startRId = startRId + blocklineinit) { long blockline = blocklineinit; if (startRId + blockline > rowcount) { @@ -345,25 +365,29 @@ ErrorCode GF3RDCreateLookTable(QString inxmlPath, QString indemPath, QString out Eigen::MatrixXd dem_y = demxyz.getData(startRId, 0, blockline, colcount, 2); Eigen::MatrixXd dem_z = demxyz.getData(startRId, 0, blockline, colcount, 3); - // 逐像素迭代计算 - double timeR = 0; + long blockrows = sar_r.rows(); long blockcols = sar_r.cols(); - double tx = 0; - double ty = 0; - double tz = 0; - double R = 0; - double slopex=0, slopey=0, slopez=0; - SatelliteOribtNode node{0,0,0,0,0,0,0,0}; - bool antflag = false; + +#pragma omp parallel for for (long i = 0; i < blockrows; i++) { + // 逐像素迭代计算 + double timeR = 0; + double tx = 0; + double ty = 0; + double tz = 0; + double R = 0; + double slopex = 0, slopey = 0, slopez = 0; + SatelliteOribtNode node{ 0,0,0,0,0,0,0,0 }; + bool antflag = false; for (long j = 0; j < blockcols; j++) { + tx = dem_x(i, j); ty = dem_y(i, j); tz = dem_z(i, j); if (RD_PSTN(refrange,lamda, timeR,R,tx,ty,tz,slopex,slopey,slopez,polyfitmodel,node,d0,d1,d2,d3,d4) == ErrorCode::SUCCESS) { - sar_r(i, j) = timeR * prf; - sar_c(i, j) = ((R - nearRange) / 2 / LIGHTSPEED) * fs; + sar_r(i, j) =( timeR+ (imagestarttime + imageendtime) / 2 -imagestarttime) * prf; + sar_c(i, j) = ((R - nearRange) * 2 / LIGHTSPEED) * fs; } else { sar_r(i, j) = -1; @@ -372,6 +396,14 @@ ErrorCode GF3RDCreateLookTable(QString inxmlPath, QString indemPath, QString out } } + // 保存结果 + omp_set_lock(&lock); + rasterRC.saveImage(sar_r, startRId, 0, 1); + rasterRC.saveImage(sar_c, startRId, 0, 2); + + omp_unset_lock(&lock); + + Eigen::MatrixXd Angle_Arr = Eigen::MatrixXd::Zero(dem_x.rows(),dem_x.cols()).array()+181; if (localincAngleFlag) { Eigen::MatrixXd demslope_x = demslope.getData(startRId, 0, blockline, colcount, 1); @@ -381,7 +413,20 @@ ErrorCode GF3RDCreateLookTable(QString inxmlPath, QString indemPath, QString out double Rst_x = 0, Rst_y = 0, Rst_z = 0, localangle = 0; double slopeR = 0; +#pragma omp parallel for for (long i = 0; i < blockrows; i++) { + // 逐像素迭代计算 + double timeR = 0; + double tx = 0; + double ty = 0; + double tz = 0; + double R = 0; + double slopex = 0, slopey = 0, slopez = 0; + SatelliteOribtNode node{ 0,0,0,0,0,0,0,0 }; + bool antflag = false; + + + for (long j = 0; j < blockcols; j++) { timeR = sar_r(i, j) / prf; slopex = demslope_x(i, j); @@ -402,16 +447,23 @@ ErrorCode GF3RDCreateLookTable(QString inxmlPath, QString indemPath, QString out Angle_Arr(i, j) = localangle; } } + + // 保存结果 + omp_set_lock(&lock); + + if (localincAngleFlag) { + localincangleRaster.saveImage(Angle_Arr, startRId, 0, 1); + } + else {} + + omp_unset_lock(&lock); + } + // 保存结果 omp_set_lock(&lock); - rasterRC.saveImage(sar_r, startRId, 0, 1); - rasterRC.saveImage(sar_c, startRId, 0, 2); - if (localincAngleFlag) { - localincangleRaster.saveImage(Angle_Arr, startRId, 0, 1); - } - else {} + processNumber = processNumber + blockrows; qDebug() << u8"\rprocess bar:\t" << processNumber * 100.0 / rowcount << " % " << "\t\t\t"; if (progressDialog.maximum() <= processNumber) { @@ -431,7 +483,8 @@ ErrorCode GF3OrthSLC( QString inRasterPath, QString inlooktablePath, QString out gdalImage slcRaster(inRasterPath);// gdalImage looktableRaster(inlooktablePath);// - gdalImage outRaster(outRasterPath);// + gdalImage outRaster = CreategdalImage(outRasterPath, looktableRaster.height, looktableRaster.width, 1, looktableRaster.gt, looktableRaster.projection, true, true);// X,Y,Z + //gdalImage outRaster(outRasterPath);// if (outRaster.height != looktableRaster.height || outRaster.width != looktableRaster.width) { qDebug() << u8"look table size is not same as outRaster size"<< looktableRaster.height <<"!="<