修正多普勒计算公式(去掉地速)
parent
d2ee22b0f1
commit
d45b397560
|
@ -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<double> 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<double> SARSimulationImageL1Dataset::getDopplerCenterCoff()
|
||||
{
|
||||
QVector<double> 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;
|
||||
}
|
||||
|
||||
|
|
|
@ -143,9 +143,11 @@ public:
|
|||
double getDopplerParametersReferenceTime();
|
||||
void setDopplerParametersReferenceTime(double v);
|
||||
|
||||
// 多普勒参数
|
||||
QVector<double> getDopplerParams();
|
||||
void setDopplerParams(double d0, double d1, double d2, double d3, double d4);
|
||||
|
||||
// 多普勒中心系数
|
||||
QVector<double> getDopplerCenterCoff();
|
||||
void setDopplerCenterCoff(double r0, double r1, double r2, double r3, double r4);
|
||||
|
||||
|
|
|
@ -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;
|
||||
}
|
||||
|
||||
|
|
|
@ -122,6 +122,7 @@
|
|||
<EnableParallelCodeGeneration>true</EnableParallelCodeGeneration>
|
||||
<LanguageStandard>stdcpp14</LanguageStandard>
|
||||
<LanguageStandard_C>stdc11</LanguageStandard_C>
|
||||
<Optimization>Disabled</Optimization>
|
||||
</ClCompile>
|
||||
</ItemDefinitionGroup>
|
||||
<ItemDefinitionGroup Condition="'$(Configuration)|$(Platform)'=='Release|ARM'">
|
||||
|
|
|
@ -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 <<"!="<<outRaster.height<<" "<<looktableRaster.width<<"!="<<outRaster.width;
|
||||
|
@ -440,7 +493,7 @@ ErrorCode GF3OrthSLC( QString inRasterPath, QString inlooktablePath, QString out
|
|||
|
||||
Eigen::MatrixXd slcImg = slcRaster.getData(0, 0, slcRaster.height, slcRaster.width, 1);
|
||||
Eigen::MatrixXd sar_r = looktableRaster.getData(0, 0, looktableRaster.height, looktableRaster.width, 1);
|
||||
Eigen::MatrixXd sar_c = looktableRaster.getData(0, 0, looktableRaster.height, looktableRaster.width, 1);
|
||||
Eigen::MatrixXd sar_c = looktableRaster.getData(0, 0, looktableRaster.height, looktableRaster.width, 2);
|
||||
Eigen::MatrixXd outImg = outRaster.getData(0, 0, outRaster.height, outRaster.width, 1);
|
||||
|
||||
long lookRows = sar_r.rows();
|
||||
|
@ -469,7 +522,7 @@ ErrorCode GF3OrthSLC( QString inRasterPath, QString inlooktablePath, QString out
|
|||
|
||||
for (long i = 0; i < lookRows; i++) {
|
||||
for (long j = 0; j < lookCols; j++) {
|
||||
p0 = {sar_r(i,j),sar_c(i,j),0};
|
||||
|
||||
lastr = std::floor(sar_r(i, j));
|
||||
nextr = std::ceil(sar_r(i, j));
|
||||
lastc = std::floor(sar_c(i, j));
|
||||
|
@ -477,11 +530,13 @@ ErrorCode GF3OrthSLC( QString inRasterPath, QString inlooktablePath, QString out
|
|||
if (lastr < 0 || lastc < 0 || nextr >= slcRows || nextc >= slcCols) {
|
||||
continue;
|
||||
}
|
||||
p11 = { sar_r(i,j),sar_c(i,j),0 };
|
||||
|
||||
p0 = { sar_c(i,j)-lastc, sar_r(i,j)-lastr,0 };
|
||||
p11 = Landpoint{ 0,0,slcImg(lastr,lastc) };
|
||||
p21 = Landpoint{ 0,1,slcImg(nextr,lastc) };
|
||||
p12 = Landpoint{ 1,0,slcImg(lastr,nextc) };
|
||||
p22 = Landpoint{ 1,1,slcImg(nextr,nextc) };
|
||||
|
||||
Bileanervalue = Bilinear_interpolation(p0, p11, p21, p12, p22);
|
||||
outImg(i, j) = Bileanervalue;
|
||||
}
|
||||
|
|
|
@ -14,10 +14,10 @@ QOrthSlrRaster::QOrthSlrRaster(QWidget *parent)
|
|||
|
||||
connect(ui.pushButtonAdd, SIGNAL(clicked(bool)), this, SLOT(onpushButtonAddClicked(bool)));
|
||||
connect(ui.pushButtonRemove, SIGNAL(clicked(bool)), this, SLOT(onpushButtonRemoveClicked(bool)));
|
||||
connect(ui.pushButtonDEMSelect, SIGNAL(clicked(bool)), this, SLOT(onpushButtonWorkSpaceClicked(bool)));
|
||||
connect(ui.pushButtonWorkSpace, SIGNAL(clicked(bool)), this, SLOT(pushButtonDEMSelectClicked(bool)));
|
||||
connect(ui.pushButtonDEMSelect, SIGNAL(clicked(bool)), this, SLOT(pushButtonDEMSelectClicked(bool)));
|
||||
connect(ui.pushButtonWorkSpace, SIGNAL(clicked(bool)), this, SLOT(onpushButtonWorkSpaceClicked(bool)));
|
||||
QObject::connect(ui.buttonBox, SIGNAL(rejected()), this, SLOT(onreject()));
|
||||
QObject::connect(ui.buttonBox, SIGNAL(accepted()), this, SLOT(onreject()));
|
||||
QObject::connect(ui.buttonBox, SIGNAL(accepted()), this, SLOT(onaccept()));
|
||||
}
|
||||
|
||||
QOrthSlrRaster::~QOrthSlrRaster()
|
||||
|
@ -61,7 +61,7 @@ void QOrthSlrRaster::onpushButtonAddClicked(bool)
|
|||
this, // 父窗口
|
||||
tr(u8"选择xml文件"), // 标题
|
||||
QString(), // 默认路径
|
||||
tr(u8"xml Files (*.xml);;All Files (*)") // Îļþ¹ýÂËÆ÷
|
||||
tr(ENVI_FILE_FORMAT_FILTER) // 文件过滤器
|
||||
);
|
||||
|
||||
// 如果用户选择了文件
|
||||
|
|
Loading…
Reference in New Issue