RasterProcessTool/GF3ProcessToolbox/GF3CalibrationAndGeocodingC...

359 lines
14 KiB
C++
Raw Normal View History

#include "GF3CalibrationAndGeocodingClass.h"
#include "SatelliteGF3xmlParser.h"
#include <QRegularExpression>
#include <QMessageBox>
ErrorCode GF3CalibrationRaster(QString inRasterPath, QString outRasterPath, double Qualifyvalue, double calibrationConst)
{
gdalImage imgraster(inRasterPath.toUtf8().constData());
if (!isExists(outRasterPath)) {
gdalImageComplex outraster = CreategdalImageComplex(outRasterPath, imgraster.height, imgraster.width, 1, imgraster.gt, imgraster.projection, false, true);
}
gdalImageComplex outraster(outRasterPath);
long blocklines = Memory1GB * 2 /8/ imgraster.width;
blocklines = blocklines < 100 ? 100 : blocklines;
Eigen::MatrixXd imgArrb1 = Eigen::MatrixXd::Zero(blocklines, imgraster.width);
Eigen::MatrixXd imgArrb2 = Eigen::MatrixXd::Zero(blocklines, imgraster.width);
Eigen::MatrixXcd imgArr = Eigen::MatrixXcd::Zero(blocklines, imgraster.width);
double quayCoff = Qualifyvalue * 1.0 / 32767;
double caliCoff = std::pow(10.0, (calibrationConst * 1.0 / 20));
qDebug() << "<EFBFBD><EFBFBD><EFBFBD><EFBFBD>ϵ<EFBFBD><EFBFBD>:\t" << quayCoff / caliCoff;
long startrow = 0;
for (startrow = 0; startrow < imgraster.height; startrow = startrow + blocklines) {
imgArrb1 = imgraster.getData(startrow, 0, blocklines, imgraster.width, 1);
imgArrb2 = imgraster.getData(startrow, 0, blocklines, imgraster.width, 2);
imgArr = outraster.getDataComplex(startrow, 0, blocklines, outraster.width, 1);
imgArr.real() = imgArrb1.array() * quayCoff / caliCoff;
imgArr.imag() = imgArrb2.array() * quayCoff / caliCoff;
outraster.saveImage(imgArr, startrow, 0, 1);
}
qDebug() << "Ӱ<EFBFBD><EFBFBD>д<EFBFBD><EFBFBD><EFBFBD>" << outRasterPath;
return ErrorCode::SUCCESS;
}
ErrorCode ImportGF3L1ARasterToWorkspace(QString inMetaxmlPath, QString inRasterPath, QString outWorkDirPath, POLARTYPEENUM polsartype)
{
SatelliteGF3xmlParser gf3xml;
gf3xml.loadFile(inMetaxmlPath);
QString filename = getFileNameWidthoutExtend(inRasterPath);
SARSimulationImageL1Dataset l1dataset;
l1dataset.OpenOrNew(outWorkDirPath, filename,gf3xml.height,gf3xml.width);
l1dataset.setSateAntPos(SatellitePos2SatelliteAntPos(gf3xml.antposs));
l1dataset.setCenterFreq(gf3xml.RadarCenterFrequency);
l1dataset.setcolCount(gf3xml.width);
l1dataset.setrowCount(gf3xml.height);
l1dataset.setNearRange(gf3xml.nearRange);
l1dataset.setRefRange(gf3xml.refRange);
l1dataset.setFs(gf3xml.eqvFs);
l1dataset.setPRF(gf3xml.eqvPRF);
l1dataset.setCenterAngle((gf3xml.incidenceAngleNearRange + gf3xml.incidenceAngleFarRange) / 2.0); // <20><><EFBFBD><EFBFBD><EFBFBD><EFBFBD>
l1dataset.setDopplerParametersReferenceTime(gf3xml.DopplerParametersReferenceTime);
l1dataset.setTotalProcessedAzimuthBandWidth(gf3xml.TotalProcessedAzimuthBandWidth);
l1dataset.setIncidenceAngleFarRange(gf3xml.incidenceAngleFarRange);
l1dataset.setIncidenceAngleNearRangeet(gf3xml.incidenceAngleNearRange);
l1dataset.setDopplerParams(gf3xml.d0, gf3xml.d1, gf3xml.d2, gf3xml.d3, gf3xml.d4);
l1dataset.setDopplerCenterCoff(gf3xml.r0, gf3xml.r1,gf3xml.r2, gf3xml.r3, gf3xml.r4);
l1dataset.setLatitudeCenter(gf3xml.latitude_center);
l1dataset.setLongitudeCenter(gf3xml.longitude_center);
l1dataset.setLatitudeTopLeft(gf3xml.latitude_topLeft);
l1dataset.setLongitudeTopLeft(gf3xml.longitude_topLeft);
l1dataset.setLatitudeTopRight(gf3xml.latitude_topRight);
l1dataset.setLongitudeTopRight(gf3xml.longitude_topRight);
l1dataset.setLatitudeBottomLeft(gf3xml.latitude_bottomLeft);
l1dataset.setLongitudeBottomLeft(gf3xml.longitude_bottomLeft);
l1dataset.setLatitudeBottomRight(gf3xml.latitude_bottomRight);
l1dataset.setLongitudeBottomRight(gf3xml.longitude_bottomRight);
l1dataset.setStartImageTime(gf3xml.start);
l1dataset.setEndImageTime(gf3xml.end);
l1dataset.saveToXml();
QString outRasterpath = l1dataset.getImageRasterPath();
ErrorCode errorcode = ErrorCode::SUCCESS;
qDebug() << "<EFBFBD><EFBFBD><EFBFBD><EFBFBD><EFBFBD><EFBFBD><EFBFBD><EFBFBD>:\t"<<inRasterPath;
switch (polsartype)
{
case POLARHH:
qDebug() << "<EFBFBD><EFBFBD><EFBFBD><EFBFBD><EFBFBD><EFBFBD>HH";
errorcode=GF3CalibrationRaster(inRasterPath, outRasterpath, gf3xml.HH_QualifyValue, gf3xml.HH_CalibrationConst);
break;
case POLARHV:
qDebug() << "<EFBFBD><EFBFBD><EFBFBD><EFBFBD><EFBFBD><EFBFBD>HH";
errorcode = GF3CalibrationRaster(inRasterPath, outRasterpath, gf3xml.HV_QualifyValue, gf3xml.HV_CalibrationConst);
break;
case POLARVH:
qDebug() << "<EFBFBD><EFBFBD><EFBFBD><EFBFBD><EFBFBD><EFBFBD>VH";
errorcode = GF3CalibrationRaster(inRasterPath, outRasterpath, gf3xml.VH_QualifyValue, gf3xml.VH_CalibrationConst);
break;
case POLARVV:
qDebug() << "<EFBFBD><EFBFBD><EFBFBD><EFBFBD><EFBFBD><EFBFBD>VV";
errorcode = GF3CalibrationRaster(inRasterPath, outRasterpath, gf3xml.VV_QualifyValue, gf3xml.VV_CalibrationConst);
break;
default:
break;
}
qDebug() << "<EFBFBD><EFBFBD><EFBFBD><EFBFBD><EFBFBD><EFBFBD><EFBFBD><EFBFBD>״̬<EFBFBD><EFBFBD>\t" << QString::fromStdString(errorCode2errInfo(errorcode));
return errorcode;
}
QVector<QString> SearchGF3DataTiff(QString inMetaxmlPath)
{
// ָ<><D6B8>·<EFBFBD><C2B7>
QString xmlpath = inMetaxmlPath; // <20>滻Ϊ<E6BBBB><CEAA><EFBFBD><EFBFBD>ʵ<EFBFBD><CAB5>·<EFBFBD><C2B7>
QString absPath = QFileInfo(xmlpath).absolutePath();
// <20><>ȡ·<C8A1><C2B7><EFBFBD><EFBFBD><EFBFBD>ڵ<EFBFBD>Ŀ¼
QDir directory(absPath);
if (!directory.exists()) {
qDebug() << "Directory does not exist:" << directory.absolutePath();
return QVector<QString>(0);
}
// <20><><EFBFBD>ù<EFBFBD><C3B9><EFBFBD><EFBFBD><EFBFBD><EFBFBD>Խ<EFBFBD><D4BD>г<EFBFBD> .tif <20>ļ<EFBFBD>
QStringList filters;
filters << "*.tif" << "*.TIF" << "*.tiff" << "*.TIFF";
QStringList files = directory.entryList(filters, QDir::Files);
qDebug() << "TIFF Files in the directory" << directory.absolutePath() << ":";
QVector<QString> filepath(0);
for (long i = 0; i < files.count(); i++) {
filepath.append(JoinPath(absPath, files[i]));
}
return filepath;
}
POLARTYPEENUM getDatasetGF3FilePolsarType(QString fileName)
{
// <20><>ȡ<EFBFBD><C8A1><EFBFBD><EFBFBD><EFBFBD><EFBFBD><EFBFBD><EFBFBD>
QRegularExpression re("_([A-Z]{2})_");
QRegularExpressionMatch match = re.match(fileName);
QString polarization = "";
if (match.hasMatch()) {
polarization = match.captured(1);
polarization = polarization.toLower().replace("_", "");
qDebug() << "Polarization extracted:" << polarization;
if (polarization == "hh") {
return POLARTYPEENUM::POLARHH;
}
else if (polarization == "hv") {
return POLARTYPEENUM::POLARHV;
}
else if (polarization == "vh") {
return POLARTYPEENUM::POLARVH;
}
else if (polarization == "vv") {
return POLARTYPEENUM::POLARVV;
}
else {
return POLARTYPEENUM::POLARUNKOWN;
}
}
else {
qDebug() << "No polarization found in the path.";
return POLARTYPEENUM::POLARUNKOWN;
}
return POLARTYPEENUM::POLARUNKOWN;
}
ErrorCode ImportGF3L1AProcess(QString inMetaxmlPath, QString outWorkDirPath)
{
QVector<QString> tiffpaths = SearchGF3DataTiff(inMetaxmlPath);
ErrorCode errorcode = ErrorCode::SUCCESS;
for (long i = 0; i < tiffpaths.count(); i++) {
QString tiffpath = tiffpaths[i];
QString filenamewithoutextern = getFileNameWidthoutExtend(tiffpath);
QString filename = getFileNameFromPath(tiffpath);
POLARTYPEENUM poltype = getDatasetGF3FilePolsarType(filename);
switch (poltype)
{
case POLARHH:
errorcode=ImportGF3L1ARasterToWorkspace(inMetaxmlPath, tiffpath, outWorkDirPath, POLARTYPEENUM::POLARHH);
break;
case POLARHV:
errorcode=ImportGF3L1ARasterToWorkspace(inMetaxmlPath, tiffpath, outWorkDirPath, POLARTYPEENUM::POLARHV);
break;
case POLARVH:
errorcode=ImportGF3L1ARasterToWorkspace(inMetaxmlPath, tiffpath, outWorkDirPath, POLARTYPEENUM::POLARVH);
break;
case POLARVV:
errorcode=ImportGF3L1ARasterToWorkspace(inMetaxmlPath, tiffpath, outWorkDirPath, POLARTYPEENUM::POLARVV);
break;
case POLARUNKOWN:
break;
default:
break;
}
if (errorcode == ErrorCode::SUCCESS) {
return errorcode;
}
else {
QMessageBox::warning(nullptr, u8"<EFBFBD><EFBFBD><EFBFBD><EFBFBD>", u8"<EFBFBD><EFBFBD><EFBFBD>ݵ<EFBFBD><EFBFBD><EFBFBD><EFBFBD><EFBFBD><EFBFBD><EFBFBD>");
break;
}
}
return errorcode;
}
ErrorCode Complex2AmpRaster(QString inComplexPath, QString outRasterPath)
{
gdalImageComplex inimg(inComplexPath);
gdalImage ampimg= CreategdalImage(outRasterPath, inimg.height, inimg.width, inimg.band_num, inimg.gt, inimg.projection, true, true);
long blocklines = Memory1GB * 2 / 8 / inimg.width;
blocklines = blocklines < 100 ? 100 : blocklines;
Eigen::MatrixXd imgArrb1 = Eigen::MatrixXd::Zero(blocklines, ampimg.width);
Eigen::MatrixXcd imgArr = Eigen::MatrixXcd::Zero(blocklines, inimg.width);
long startrow = 0;
for (startrow = 0; startrow < inimg.height; startrow = startrow + blocklines) {
imgArrb1 = ampimg.getData(startrow, 0, blocklines, inimg.width, 1);
imgArr = inimg.getData(startrow, 0, blocklines, inimg.width, 2);
imgArrb1 = imgArr.array().abs();
ampimg.saveImage(imgArrb1, startrow, 0, 1);
}
qDebug() << "Ӱ<EFBFBD><EFBFBD>д<EFBFBD><EFBFBD><EFBFBD>" << outRasterPath;
return ErrorCode::SUCCESS;
}
ErrorCode Complex2PhaseRaster(QString inComplexPath, QString outRasterPath)
{
gdalImageComplex inimg(inComplexPath);
gdalImage ampimg = CreategdalImage(outRasterPath, inimg.height, inimg.width, inimg.band_num, inimg.gt, inimg.projection, true, true);
long blocklines = Memory1GB * 2 / 8 / inimg.width;
blocklines = blocklines < 100 ? 100 : blocklines;
Eigen::MatrixXd imgArrb1 = Eigen::MatrixXd::Zero(blocklines, ampimg.width);
Eigen::MatrixXcd imgArr = Eigen::MatrixXcd::Zero(blocklines, inimg.width);
long startrow = 0;
for (startrow = 0; startrow < inimg.height; startrow = startrow + blocklines) {
imgArrb1 = ampimg.getData(startrow, 0, blocklines, inimg.width, 1);
imgArr = inimg.getData(startrow, 0, blocklines, inimg.width, 2);
imgArrb1 = imgArr.array().arg();
ampimg.saveImage(imgArrb1, startrow, 0, 1);
}
qDebug() << "Ӱ<EFBFBD><EFBFBD>д<EFBFBD><EFBFBD><EFBFBD>" << outRasterPath;
return ErrorCode::SUCCESS;
}
ErrorCode Complex2dBRaster(QString inComplexPath, QString outRasterPath)
{
gdalImageComplex inimg(inComplexPath);
gdalImage ampimg = CreategdalImage(outRasterPath, inimg.height, inimg.width, inimg.band_num, inimg.gt, inimg.projection, true, true);
long blocklines = Memory1GB * 2 / 8 / inimg.width;
blocklines = blocklines < 100 ? 100 : blocklines;
Eigen::MatrixXd imgArrb1 = Eigen::MatrixXd::Zero(blocklines, ampimg.width);
Eigen::MatrixXcd imgArr = Eigen::MatrixXcd::Zero(blocklines, inimg.width);
long startrow = 0;
for (startrow = 0; startrow < inimg.height; startrow = startrow + blocklines) {
imgArrb1 = ampimg.getData(startrow, 0, blocklines, inimg.width, 1);
imgArr = inimg.getData(startrow, 0, blocklines, inimg.width, 2);
imgArrb1 = imgArr.array().abs().log10()*20.0;
ampimg.saveImage(imgArrb1, startrow, 0, 1);
}
qDebug() << "Ӱ<EFBFBD><EFBFBD>д<EFBFBD><EFBFBD><EFBFBD>" << outRasterPath;
return ErrorCode::SUCCESS;
}
ErrorCode ResampleDEM(QString indemPath, QString outdemPath,double gridx,double gridy)
{
double gridlat = gridy / 110000.0;
double gridlon = gridx / 100000.0;
long espgcode = GetEPSGFromRasterFile(indemPath.toUtf8().constData());
if (espgcode == 4326) {
resampleRaster(indemPath.toUtf8().constData(), outdemPath.toUtf8().constData(), gridlon, gridlat);
return ErrorCode::SUCCESS;
}
else {
QMessageBox::warning(nullptr,u8"<EFBFBD><EFBFBD><EFBFBD><EFBFBD>",u8"<EFBFBD><EFBFBD><EFBFBD><EFBFBD><EFBFBD><EFBFBD>WGS84<EFBFBD><EFBFBD><EFBFBD><EFBFBD><EFBFBD><EFBFBD>DEMӰ<EFBFBD><EFBFBD>");
return ErrorCode::FAIL;
}
}
ErrorCode GF3RDCreateLookTable(QString inxmlPath, QString indemPath,QString outworkdir, QString outlooktablePath)
{
QString indemfilename = getFileNameWidthoutExtend(indemPath);
QString resampleDEMXYZPath = JoinPath(outworkdir, indemfilename + "_XYZ.tif");
QString resampleDEMSLOPEPath = JoinPath(outworkdir, indemfilename + "_slopevector.tif");
DEM2XYZRasterAndSlopRaster(indemPath, resampleDEMXYZPath, resampleDEMSLOPEPath);
gdalImage demxyz(resampleDEMXYZPath);// X,Y,Z
gdalImage demslope(resampleDEMSLOPEPath);// X,Y,Z
gdalImage rasterRC = CreategdalImage(outlooktablePath, demxyz.height, demxyz.width,2, demxyz.gt, demxyz.projection, true, true);// X,Y,Z
SARSimulationImageL1Dataset l1dataset;
l1dataset.Open(inxmlPath);
// RD <20><><EFBFBD>ز<EFBFBD><D8B2><EFBFBD>
QVector<double> dopplers = l1dataset.getDopplerParams();
double d0 = dopplers[0];
double d1 = dopplers[1];
double d2 = dopplers[2];
double d3 = dopplers[3];
double d4 = dopplers[4];
double fs = l1dataset.getFs();
double prf = (l1dataset.getEndImageTime()-l1dataset.getStartImageTime())/(l1dataset.getrowCount()-1);
// <20><><EFBFBD><EFBFBD><EFBFBD><EFBFBD><EFBFBD><EFBFBD>ģ<EFBFBD><C4A3>
return ErrorCode::SUCCESS;
}
ErrorCode GF3RDCreateLookTable(QString inxmlPath, QString indemPath, QString outworkdir, double gridx, double gridy)
{
SARSimulationImageL1Dataset l1dataset;
l1dataset.Open(inxmlPath);
DemBox box = l1dataset.getExtend();
double dlon = 0.1*(box.max_lon - box.min_lon);
double dlat = 0.1 * (box.max_lat - box.min_lat);
double minlat = box.min_lat - dlat;
double minlon = box.min_lon - dlon;
double maxlat = box.max_lat + dlat;
double maxlon = box.max_lon + dlon;
// <20>ü<EFBFBD>DEM
QString filename = getFileNameWidthoutExtend(l1dataset.getxmlFilePath());
filename = filename.right(5);
QString clipDEMPath = JoinPath(outworkdir, getFileNameWidthoutExtend(indemPath) + filename+"_clip.tif");
cropRasterByLatLon(indemPath.toUtf8().constData(), clipDEMPath.toUtf8().constData(), minlon, maxlon, minlat, maxlat);
QString resampleDEMPath = JoinPath(outworkdir, getFileNameWidthoutExtend(indemPath) + filename + "_clip_resample.tif");
resampleRaster(clipDEMPath.toUtf8().constData(), resampleDEMPath.toUtf8().constData(), gridx, gridy);
QString outlooktablePath= JoinPath(outworkdir, getFileNameWidthoutExtend(l1dataset.getxmlFilePath()) + "_looktable.tif");
if (GF3RDCreateLookTable(inxmlPath, resampleDEMPath, outworkdir, outlooktablePath) != ErrorCode::SUCCESS) {
qDebug() << "<EFBFBD><EFBFBD><EFBFBD>ұ<EFBFBD><EFBFBD><EFBFBD><EFBFBD>ɴ<EFBFBD><EFBFBD><EFBFBD><EFBFBD><EFBFBD>\t"+ getFileNameWidthoutExtend(inxmlPath);
return ErrorCode::FAIL;
}
return ErrorCode::SUCCESS;
}