RasterProcessTool/SimulationSAR/RTPCProcessCls.cpp

1100 lines
43 KiB
C++
Raw Blame History

This file contains ambiguous Unicode characters!

This file contains ambiguous Unicode characters that may be confused with others in your current locale. If your use case is intentional and legitimate, you can safely ignore this warning. Use the Escape button to highlight these characters.

#include "stdafx.h"
#include "RTPCProcessCls.h"
#include "BaseConstVariable.h"
#include "SARSatelliteSimulationAbstractCls.h"
#include "SARSimulationTaskSetting.h"
#include "SatelliteOribtModel.h"
#include <QDebug>
#include "ImageOperatorBase.h"
#include "GeoOperator.h"
#include "EchoDataFormat.h"
#include <QDir>
#include <QDatetime>
#include <omp.h>
#include <QProgressDialog>
#include <QMessageBox>
#ifdef DEBUGSHOWDIALOG
#include "ImageShowDialogClass.h"
#endif
#ifdef __CUDANVCC___
#include "GPUTool.cuh"
#endif // __CUDANVCC___
#include <Imageshow/ImageShowDialogClass.h>
RTPCProcessCls::RTPCProcessCls()
{
this->PluseCount = 0;
this->PlusePoint = 0;
this->TaskSetting = nullptr;
this->EchoSimulationData = nullptr;
this->DEMTiffPath = "";
this->LandCoverPath = "";
this->HHSigmaPath = "";
this->HVSigmaPath = "";
this->VHSigmaPath = "";
this->VVSigmaPath = "";
this->OutEchoPath = "";
this->DEMTiffPath.clear();
this->LandCoverPath.clear();
this->HHSigmaPath.clear();
this->HVSigmaPath.clear();
this->VHSigmaPath.clear();
this->VVSigmaPath.clear();
this->OutEchoPath.clear();
this->SigmaDatabasePtr = std::shared_ptr<SigmaDatabase>(new SigmaDatabase);
}
RTPCProcessCls::~RTPCProcessCls()
{
}
void RTPCProcessCls::setTaskSetting(std::shared_ptr < AbstractSARSatelliteModel> TaskSetting)
{
this->TaskSetting = std::shared_ptr < AbstractSARSatelliteModel>(TaskSetting);
qDebug() << "RTPCProcessCls::setTaskSetting";
}
void RTPCProcessCls::setEchoSimulationDataSetting(std::shared_ptr<EchoL0Dataset> EchoSimulationData)
{
this->EchoSimulationData = std::shared_ptr<EchoL0Dataset>(EchoSimulationData);
qDebug() << "RTPCProcessCls::setEchoSimulationDataSetting";
}
void RTPCProcessCls::setTaskFileName(QString EchoFileName)
{
this->TaskFileName = EchoFileName;
}
void RTPCProcessCls::setDEMTiffPath(QString DEMTiffPath)
{
this->DEMTiffPath = DEMTiffPath;
}
void RTPCProcessCls::setLandCoverPath(QString LandCoverPath)
{
this->LandCoverPath = LandCoverPath;
}
void RTPCProcessCls::setHHSigmaPath(QString HHSigmaPath)
{
this->HHSigmaPath = HHSigmaPath;
}
void RTPCProcessCls::setHVSigmaPath(QString HVSigmaPath)
{
this->HVSigmaPath = HVSigmaPath;
}
void RTPCProcessCls::setVHSigmaPath(QString VHSigmaPath)
{
this->VHSigmaPath = VHSigmaPath;
}
void RTPCProcessCls::setVVSigmaPath(QString VVSigmaPath)
{
this->VVSigmaPath = VVSigmaPath;
}
void RTPCProcessCls::setOutEchoPath(QString OutEchoPath)
{
this->OutEchoPath = OutEchoPath;
}
ErrorCode RTPCProcessCls::Process(long num_thread)
{
// RTPC <20>
qDebug() << u8"params init ....";
ErrorCode stateCode = this->InitParams();
if (stateCode != ErrorCode::SUCCESS) {
return stateCode;
}
else {}
qDebug() << "DEMMainProcess";
stateCode = this->DEMPreprocess();
if (stateCode != ErrorCode::SUCCESS) {
return stateCode;
}
else {}
qDebug() << "RTPCMainProcess";
//stateCode = this->RTPCMainProcess(num_thread);
stateCode = this->RTPCMainProcess_GPU( );
if (stateCode != ErrorCode::SUCCESS) {
return stateCode;
}
else {}
return ErrorCode::SUCCESS;
}
ErrorCode RTPCProcessCls::InitParams()
{
if (nullptr == this->TaskSetting || this->DEMTiffPath.isEmpty() ||
this->LandCoverPath.isEmpty() || this->HHSigmaPath.isEmpty() ||
this->HVSigmaPath.isEmpty() || this->VHSigmaPath.isEmpty() ||
this->VVSigmaPath.isEmpty()) {
return ErrorCode::RTPC_PARAMSISEMPTY;
}
else {
}
// <20><>һ<EFBFBD><D2BB><EFBFBD><EFBFBD><EFBFBD><EFBFBD>·<EFBFBD><C2B7>
this->OutEchoPath = QDir(this->OutEchoPath).absolutePath();
// <20>ز<EFBFBD><D8B2><EFBFBD>С
double imgStart_end = this->TaskSetting->getSARImageEndTime() - this->TaskSetting->getSARImageStartTime();
this->PluseCount = ceil(imgStart_end * this->TaskSetting->getPRF());
double rangeTimeSample = (this->TaskSetting->getFarRange() - this->TaskSetting->getNearRange()) * 2.0 / LIGHTSPEED;
this->PlusePoint = ceil(rangeTimeSample * this->TaskSetting->getFs());
// <20><>ʼ<EFBFBD><CABC><EFBFBD>ز<EFBFBD><D8B2><EFBFBD><EFBFBD><EFBFBD>λ<EFBFBD><CEBB>
qDebug() << "--------------Echo Data Setting ---------------------------------------";
this->EchoSimulationData = std::shared_ptr<EchoL0Dataset>(new EchoL0Dataset);
this->EchoSimulationData->setCenterFreq(this->TaskSetting->getCenterFreq());
this->EchoSimulationData->setNearRange(this->TaskSetting->getNearRange());
this->EchoSimulationData->setFarRange(this->TaskSetting->getFarRange());
this->EchoSimulationData->setFs(this->TaskSetting->getFs());
this->EchoSimulationData->setCenterAngle(this->TaskSetting->getCenterLookAngle());
this->EchoSimulationData->setLookSide(this->TaskSetting->getIsRightLook() ? "R" : "L");
this->EchoSimulationData->OpenOrNew(OutEchoPath, TaskFileName, PluseCount, PlusePoint);
QString tmpfolderPath = QDir(OutEchoPath).filePath("tmp");
if (QDir(tmpfolderPath).exists() == false) {
QDir(OutEchoPath).mkpath(tmpfolderPath);
}
this->tmpfolderPath = tmpfolderPath;
return ErrorCode::SUCCESS;
}
ErrorCode RTPCProcessCls::DEMPreprocess()
{
this->demxyzPath = QDir(tmpfolderPath).filePath("demxyz.tif");
gdalImage demds(this->DEMTiffPath);
gdalImage demxyz = CreategdalImage(demxyzPath, demds.height, demds.width, 3, demds.gt, demds.projection, true, true);// X,Y,Z
// <20>ֿ<EFBFBD><D6BF><EFBFBD><EFBFBD>㲢ת<E3B2A2><D7AA>Ϊ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 };
double R = 0;
double dem_row = 0, dem_col = 0, dem_alt = 0;
long line_invert = 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;
Eigen::MatrixXd xyzdata_z = demdata.array() * 0;
int datarows = demdata.rows();
int datacols = demdata.cols();
for (int i = 0; i < datarows; i++) {
for (int j = 0; j < datacols; j++) {
rowidx = i + max_rows_ids;
colidx = j;
demds.getLandPoint(rowidx, colidx, demdata(i, j), LandP); // <20><>ȡ<EFBFBD><C8A1><EFBFBD><EFBFBD><EFBFBD><EFBFBD><EFBFBD><EFBFBD>
LLA2XYZ(LandP, GERpoint); // <20><>γ<EFBFBD><CEB3>ת<EFBFBD><D7AA>Ϊ<EFBFBD><CEAA><EFBFBD><EFBFBD><EFBFBD><EFBFBD><EFBFBD><EFBFBD>ϵ
xyzdata_x(i, j) = GERpoint.x;
xyzdata_y(i, j) = GERpoint.y;
xyzdata_z(i, j) = GERpoint.z;
}
}
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);
}
// <20><><EFBFBD><EFBFBD><EFBFBD><EFBFBD><EFBFBD><EFBFBD><EFBFBD><EFBFBD>
this->demsloperPath = QDir(tmpfolderPath).filePath("demsloper.tif");
this->demmaskPath = QDir(tmpfolderPath).filePath("demmask.tif");
gdalImage demsloperxyz = CreategdalImage(this->demsloperPath, demds.height, demds.width, 4, demds.gt, demds.projection, true, true);// X,Y,Z,cosangle
gdalImage demmask = CreategdalImage(this->demmaskPath, demxyz.height, demxyz.width, 1, demxyz.gt, demxyz.projection, true, true);// X,Y,Z
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 maskdata = demmask.getData(start_ids - 1, 0, line_invert + 2, demxyz.width, 1);
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);
Eigen::MatrixXd demsloper_angle = demsloperxyz.getData(start_ids - 1, 0, line_invert + 2, demxyz.width, 4);
maskdata = maskdata.array() * 0;
Landpoint p0, p1, p2, p3, p4, pslopeVector, pp;
Vector3D slopeVector;
dem_rows = maskdata.rows();
dem_cols = maskdata.cols();
double sloperAngle = 0;
Vector3D Zaxis = { 0,0,1 };
double rowidx = 0, colidx = 0;
for (long i = 1; i < dem_rows - 1; i++) {
for (long j = 1; j < dem_cols - 1; j++) {
rowidx = i + startlineid;
colidx = j;
demds.getLandPoint(rowidx, colidx, demdata(i, j), p0);
demds.getLandPoint(rowidx - 1, colidx, demdata(i - 1, j), p1);
demds.getLandPoint(rowidx, colidx - 1, demdata(i, j - 1), p2);
demds.getLandPoint(rowidx + 1, colidx, demdata(i + 1, j), p3);
demds.getLandPoint(rowidx, colidx + 1, demdata(i, j + 1), p4);
pslopeVector = getSlopeVector(p0, p1, p2, p3, p4); // <20><><EFBFBD><EFBFBD><EFBFBD><EFBFBD><EFBFBD><EFBFBD>ʸ<EFBFBD><CAB8>
slopeVector = { pslopeVector.lon,pslopeVector.lat,pslopeVector.ati };
pp = LLA2XYZ(p0);
Zaxis.x = pp.lon;
Zaxis.y = pp.lat;
Zaxis.z = pp.ati;
sloperAngle = getCosAngle(slopeVector, Zaxis); // <20><><EFBFBD><EFBFBD><EFBFBD><EFBFBD><EFBFBD><EFBFBD><EFBFBD><EFBFBD>
demsloper_x(i, j) = slopeVector.x;
demsloper_y(i, j) = slopeVector.y;
demsloper_z(i, j) = slopeVector.z;
demsloper_angle(i, j) = sloperAngle;
maskdata(i, j)++;
}
}
demmask.saveImage(maskdata, start_ids - 1, 0, 1);
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;
}
ErrorCode RTPCProcessCls::RTPCMainProcess_GPU( )
{
double widthSpace = LIGHTSPEED / 2 / this->TaskSetting->getFs();
double prf_time = 0;
double dt = 1 / this->TaskSetting->getPRF();// <20><>ȡÿ<C8A1><C3BF><EFBFBD><EFBFBD><EFBFBD><EFBFBD><EFBFBD><EFBFBD>ʱ<EFBFBD><CAB1><EFBFBD><EFBFBD><EFBFBD><EFBFBD>
bool antflag = true; // <20><><EFBFBD><EFBFBD><EFBFBD><EFBFBD><EFBFBD>߷<EFBFBD><DFB7><EFBFBD>ͼ
Landpoint LandP{ 0,0,0 };
Point3 GERpoint{ 0,0,0 };
double R = 0;
double dem_row = 0, dem_col = 0, dem_alt = 0;
long double imageStarttime = 0;
imageStarttime = this->TaskSetting->getSARImageStartTime();
//std::vector<SatelliteOribtNode> sateOirbtNodes(this->PluseCount);
std::shared_ptr<SatelliteOribtNode[]> sateOirbtNodes(new SatelliteOribtNode[this->PluseCount], delArrPtr);
{ // <20><>̬<EFBFBD><CCAC><EFBFBD>㲻ͬ
qDebug() << "Ant position finished started !!!";
// <20><><EFBFBD><EFBFBD><EFBFBD><EFBFBD>̬
std::shared_ptr<double> antpos = this->EchoSimulationData->getAntPos();
double dAt = 1e-6;
double prf_time_dt = 0;
Landpoint InP{ 0,0,0 }, outP{ 0,0,0 };
for (long prf_id = 0; prf_id < this->PluseCount; prf_id++) {
prf_time = dt * prf_id;
prf_time_dt = prf_time + dAt;
SatelliteOribtNode sateOirbtNode;
SatelliteOribtNode sateOirbtNode_dAt;
this->TaskSetting->getSatelliteOribtNode(prf_time, sateOirbtNode, antflag);
this->TaskSetting->getSatelliteOribtNode(prf_time_dt, sateOirbtNode_dAt, antflag);
sateOirbtNode.AVx = (sateOirbtNode_dAt.Vx - sateOirbtNode.Vx) / dAt; // <20><><EFBFBD>ٶ<EFBFBD>
sateOirbtNode.AVy = (sateOirbtNode_dAt.Vy - sateOirbtNode.Vy) / dAt;
sateOirbtNode.AVz = (sateOirbtNode_dAt.Vz - sateOirbtNode.Vz) / dAt;
InP.lon = sateOirbtNode.Px;
InP.lat = sateOirbtNode.Py;
InP.ati = sateOirbtNode.Pz;
outP = XYZ2LLA(InP);
antpos.get()[prf_id * 19 + 0] = prf_time + imageStarttime;
antpos.get()[prf_id * 19 + 1] = sateOirbtNode.Px;
antpos.get()[prf_id * 19 + 2] = sateOirbtNode.Py;
antpos.get()[prf_id * 19 + 3] = sateOirbtNode.Pz;
antpos.get()[prf_id * 19 + 4] = sateOirbtNode.Vx;
antpos.get()[prf_id * 19 + 5] = sateOirbtNode.Vy;
antpos.get()[prf_id * 19 + 6] = sateOirbtNode.Vz;
antpos.get()[prf_id * 19 + 7] = sateOirbtNode.AntDirecX;
antpos.get()[prf_id * 19 + 8] = sateOirbtNode.AntDirecY;
antpos.get()[prf_id * 19 + 9] = sateOirbtNode.AntDirecZ;
antpos.get()[prf_id * 19 + 10] = sateOirbtNode.AVx;
antpos.get()[prf_id * 19 + 11] = sateOirbtNode.AVy;
antpos.get()[prf_id * 19 + 12] = sateOirbtNode.AVz;
antpos.get()[prf_id * 19 + 13] = sateOirbtNode.zeroDopplerDirectX;
antpos.get()[prf_id * 19 + 14] = sateOirbtNode.zeroDopplerDirectY;
antpos.get()[prf_id * 19 + 15] = sateOirbtNode.zeroDopplerDirectZ;
antpos.get()[prf_id * 19 + 16] = outP.lon;
antpos.get()[prf_id * 19 + 17] = outP.lat;
antpos.get()[prf_id * 19 + 18] = outP.ati;
sateOirbtNodes[prf_id] = sateOirbtNode;
}
this->EchoSimulationData->saveAntPos(antpos);
antpos.reset();
qDebug() << "Ant position finished sucessfully !!!";
}
// <20>ز<EFBFBD>
long echoIdx = 0;
double NearRange = this->EchoSimulationData->getNearRange(); // <20><>б<EFBFBD><D0B1>
double FarRange = this->EchoSimulationData->getFarRange();
double TimgNearRange = 2 * NearRange / LIGHTSPEED;
double TimgFarRange = 2 * FarRange / LIGHTSPEED;
double Fs = this->TaskSetting->getFs(); // <20><><EFBFBD><EFBFBD><EFBFBD><EFBFBD><EFBFBD><EFBFBD><EFBFBD><EFBFBD><EFBFBD><EFBFBD>
double Pt = this->TaskSetting->getPt() * this->TaskSetting->getGri();// <20><><EFBFBD><EFBFBD><EFBFBD><EFBFBD>ѹ 1v
//double GainAntLen = -3;// -3dB Ϊ<><CEAA><EFBFBD>߰뾶
long pluseCount = this->PluseCount;
double lamda = this->TaskSetting->getCenterLamda(); // <20><><EFBFBD><EFBFBD>
// <20><><EFBFBD>߷<EFBFBD><DFB7><EFBFBD>ͼ
std::shared_ptr<AbstractRadiationPattern> TransformPattern = this->TaskSetting->getTransformRadiationPattern(); // <20><><EFBFBD><EFBFBD><EFBFBD><EFBFBD><EFBFBD>߷<EFBFBD><DFB7><EFBFBD>ͼ
std::shared_ptr<AbstractRadiationPattern> ReceivePattern = this->TaskSetting->getReceiveRadiationPattern(); // <20><><EFBFBD><EFBFBD><EFBFBD><EFBFBD><EFBFBD>߷<EFBFBD><DFB7><EFBFBD>ͼ
long PlusePoint = this->EchoSimulationData->getPlusePoints();
// <20><>ʼ<EFBFBD><CABC><EFBFBD>ز<EFBFBD>
this->EchoSimulationData->initEchoArr(std::complex<double>(0, 0));
POLARTYPEENUM polartype = this->TaskSetting->getPolarType();
#ifndef __CUDANVCC___
QMessageBox::information(this, u8"<EFBFBD><EFBFBD><EFBFBD><EFBFBD><EFBFBD><EFBFBD>ʾ", u8"<EFBFBD><EFBFBD>ȷ<EFBFBD><EFBFBD><EFBFBD><EFBFBD>װ<EFBFBD><EFBFBD>CUDA<EFBFBD><EFBFBD>");
#else
// RTPC CUDA<44>
if (pluseCount * 4 * 18 > Memory1MB * 100) {
long max = Memory1MB * 100 / 4 / 20 / PluseCount;
QMessageBox::warning(nullptr, u8"<EFBFBD><EFBFBD><EFBFBD><EFBFBD><EFBFBD>̫<EFBFBD><EFBFBD><EFBFBD><EFBFBD>", u8"<EFBFBD><EFBFBD>ǰƵ<EFBFBD><EFBFBD><EFBFBD><EFBFBD><EFBFBD>£<EFBFBD><EFBFBD><EFBFBD><EFBFBD><EFBFBD><EFBFBD><EFBFBD><EFBFBD><EFBFBD><EFBFBD><EFBFBD><EFBFBD><EFBFBD>Ϊ<EFBFBD><EFBFBD>" + QString::number(max));
}
gdalImage demxyz(this->demxyzPath);// <20><><EFBFBD><EFBFBD><EFBFBD><EFBFBD><EFBFBD><EFBFBD><EFBFBD><EFBFBD>
gdalImage demlandcls(this->LandCoverPath);// <20>ر<EFBFBD><D8B1><EFBFBD><EFBFBD><EFBFBD><EFBFBD><EFBFBD><EFBFBD><EFBFBD>
gdalImage demsloperxyz(this->demsloperPath);// <20><><EFBFBD><EFBFBD><EFBFBD><EFBFBD><EFBFBD><EFBFBD>
// <20><><EFBFBD><EFBFBD><EFBFBD><EFBFBD><EFBFBD>ֿ<EFBFBD><D6BF><EFBFBD><EFBFBD><EFBFBD>
long demRow = demxyz.height;
long demCol = demxyz.width;
long blokline = 100;
// ÿ<><C3BF> 250MB*16 = 4GB
blokline = Memory1MB * 500 / 8 / demCol;
blokline = blokline < 1 ? 1 : blokline;
bool bloklineflag = false;
// <20><><EFBFBD><EFBFBD><EFBFBD><EFBFBD><EFBFBD><EFBFBD><EFBFBD><EFBFBD><EFBFBD>߷<EFBFBD><DFB7><EFBFBD>ͼ
double Tminphi = TransformPattern->getMinPhi();
double Tmaxphi = TransformPattern->getMaxPhi();
double Tmintheta = TransformPattern->getMinTheta();
double Tmaxtheta = TransformPattern->getMaxTheta();
long Tphinum = TransformPattern->getPhis().size();
long Tthetanum = TransformPattern->getThetas().size();
double TstartTheta = Tmintheta;
double TstartPhi = Tminphi;
double Tdtheta = (Tmaxtheta - Tmintheta) / (Tthetanum - 1);
double Tdphi = (Tmaxphi - Tminphi) / (Tphinum - 1);
float* h_TantPattern = (float*)mallocCUDAHost(sizeof(float) * Tthetanum * Tphinum);
float* d_TantPattern = (float*)mallocCUDADevice(sizeof(float) * Tthetanum * Tphinum);
for (long i = 0; i < Tthetanum; i++) {
for (long j = Tphinum - 1; j >=0 ; j--) {
//h_TantPattern[i * Tphinum + j] = TransformPattern->getGainLearThetaPhi(TstartTheta + i * Tdtheta, TstartPhi + j * Tdphi);
h_TantPattern[i * Tphinum + j] = TransformPattern->getGain(TstartTheta + i * Tdtheta, TstartPhi + j * Tdphi);
}
}
testOutAntPatternTrans("TransPattern.bin", h_TantPattern, TstartTheta, Tdtheta, TstartPhi, Tdphi, Tthetanum, Tphinum);
for (long i = 0; i < Tthetanum; i++) {
for (long j = 0; j < Tphinum; j++) {
h_TantPattern[i * Tphinum + j] = powf(10.0, h_TantPattern[i * Tphinum + j]/10);
}
}
HostToDevice(h_TantPattern, d_TantPattern, sizeof(float)* Tthetanum* Tphinum);
// <20><><EFBFBD><EFBFBD><EFBFBD><EFBFBD><EFBFBD><EFBFBD><EFBFBD><EFBFBD><EFBFBD>߷<EFBFBD><DFB7><EFBFBD>ͼ
double Rminphi = ReceivePattern->getMinPhi();
double Rmaxphi = ReceivePattern->getMaxPhi();
double Rmintheta = ReceivePattern->getMinTheta();
double Rmaxtheta = ReceivePattern->getMaxTheta();
long Rphinum = ReceivePattern->getPhis().size();
long Rthetanum = ReceivePattern->getThetas().size();
double RstartTheta = Rmintheta;
double RstartPhi = Rminphi;
double Rdtheta = (Rmaxtheta - Rmintheta) / (Rthetanum - 1);
double Rdphi = (Rmaxphi - Rminphi) / (Rphinum - 1);
float* h_RantPattern = (float*)mallocCUDAHost(sizeof(float) * Rthetanum * Rphinum);
float* d_RantPattern = (float*)mallocCUDADevice(sizeof(float) * Rthetanum * Rphinum);
for (long i = 0; i < Rthetanum; i++) {
for (long j = 0; j < Rphinum; j++) {
//h_RantPattern[i * Rphinum + j] = ReceivePattern->getGainLearThetaPhi(RstartTheta + i * Rdtheta, RstartPhi + j * Rdphi);
h_RantPattern[i * Rphinum + j] = ReceivePattern->getGain(RstartTheta + i * Rdtheta, RstartPhi + j * Rdphi);
}
}
testOutAntPatternTrans("ReceivePattern.bin", h_RantPattern, Rmintheta, Rdtheta, RstartPhi, Rdphi, Rthetanum, Rphinum);
for (long i = 0; i < Tthetanum; i++) {
for (long j = 0; j < Tphinum; j++) {
h_RantPattern[i * Tphinum + j] = powf(10.0, h_RantPattern[i * Tphinum + j] / 10);
}
}
HostToDevice(h_RantPattern, d_RantPattern, sizeof(float) * Rthetanum * Rphinum);
//<2F><><EFBFBD><EFBFBD><EFBFBD>ر<EFBFBD><D8B1><EFBFBD><EFBFBD><EFBFBD>
QMap<long, long> clamap;
long clamapid = 0;
long startline = 0;
for (startline = 0; startline < demRow; startline = startline + blokline) {
Eigen::MatrixXd clsland = demlandcls.getData(startline, 0, blokline, demlandcls.width, 1);
long clsrows = clsland.rows();
long clscols = clsland.cols();
long clsid = 0;
for (long ii = 0; ii < clsrows; ii++) {
for (long jj = 0; jj < clscols; jj++) {
clsid = clsland(ii, jj);
if (clamap.contains(clsid)) {}
else {
clamap.insert(clsid, clamapid);
clamapid = clamapid + 1;
}
}
}
}
std::cout << "class id recoding" << std::endl;
for (long id : clamap.keys()) {
std::cout << id << " -> " << clamap[id] << std::endl;
}
CUDASigmaParam* h_clsSigmaParam = (CUDASigmaParam*)mallocCUDAHost(sizeof(CUDASigmaParam) * clamapid);
CUDASigmaParam* d_clsSigmaParam = (CUDASigmaParam*)mallocCUDADevice(sizeof(CUDASigmaParam) * clamapid);
{
std::map<long, SigmaParam> tempSigmaParam = this->SigmaDatabasePtr->getsigmaParams(polartype);
for (long id : clamap.keys()) {
SigmaParam tempp = tempSigmaParam[id];
h_clsSigmaParam[clamap[id]].p1 = tempp.p1;
h_clsSigmaParam[clamap[id]].p2 = tempp.p2;
h_clsSigmaParam[clamap[id]].p3 = tempp.p3;
h_clsSigmaParam[clamap[id]].p4 = tempp.p4;
h_clsSigmaParam[clamap[id]].p5 = tempp.p5;
h_clsSigmaParam[clamap[id]].p6 = tempp.p6;
}
// <20><>ӡ<EFBFBD><D3A1>־
std::cout << "sigma params:" << std::endl;
std::cout << "classid:\tp1\tp2\tp3\tp4\tp5\tp6"<<std::endl;
for (long ii = 0; ii < clamapid; ii++) {
std::cout << ii << ":\t" << h_clsSigmaParam[ii].p1;
std::cout << "\t" << h_clsSigmaParam[ii].p2;
std::cout << "\t" << h_clsSigmaParam[ii].p3;
std::cout << "\t" << h_clsSigmaParam[ii].p4;
std::cout << "\t" << h_clsSigmaParam[ii].p5;
std::cout << "\t" << h_clsSigmaParam[ii].p6<<std::endl;
}
std::cout << "";
}
HostToDevice(h_clsSigmaParam, d_clsSigmaParam, sizeof(CUDASigmaParam) * clamapid);
Eigen::MatrixXd dem_x = demxyz.getData(0, 0, blokline, demxyz.width, 1); // <20><><EFBFBD><EFBFBD><EFBFBD><EFBFBD><EFBFBD><EFBFBD>
long tempDemRows = dem_x.rows();
long tempDemCols = dem_x.cols();
Eigen::MatrixXd dem_y = Eigen::MatrixXd::Zero(tempDemRows, tempDemCols);
Eigen::MatrixXd dem_z = Eigen::MatrixXd::Zero(tempDemRows, tempDemCols);
Eigen::MatrixXd demsloper_x = Eigen::MatrixXd::Zero(tempDemRows, tempDemCols);
Eigen::MatrixXd demsloper_y = Eigen::MatrixXd::Zero(tempDemRows, tempDemCols);
Eigen::MatrixXd demsloper_z = Eigen::MatrixXd::Zero(tempDemRows, tempDemCols);
Eigen::MatrixXd sloperAngle = Eigen::MatrixXd::Zero(tempDemRows, tempDemCols);
float* h_dem_x;
float* h_dem_y;
float* h_dem_z;
float* h_demsloper_x;
float* h_demsloper_y;
float* h_demsloper_z;
float* h_demsloper_angle;
float* d_dem_x;
float* d_dem_y;
float* d_dem_z;
float* d_demsloper_x;
float* d_demsloper_y;
float* d_demsloper_z;
float* d_demsloper_angle;
h_dem_x=(float* )mallocCUDAHost( sizeof(float) * blokline * tempDemCols);
h_dem_y=(float* )mallocCUDAHost(sizeof(float) * blokline * tempDemCols);
h_dem_z=(float* )mallocCUDAHost( sizeof(float) * blokline * tempDemCols);
h_demsloper_x=(float* )mallocCUDAHost( sizeof(float) * blokline * tempDemCols);
h_demsloper_y=(float* )mallocCUDAHost( sizeof(float) * blokline * tempDemCols);
h_demsloper_z=(float* )mallocCUDAHost(sizeof(float) * blokline * tempDemCols);
h_demsloper_angle= (float*)mallocCUDAHost( sizeof(float) * blokline * tempDemCols);
d_dem_x=(float* )mallocCUDADevice( sizeof(float) * blokline * tempDemCols); // 7
d_dem_y=(float* )mallocCUDADevice( sizeof(float) * blokline * tempDemCols);
d_dem_z=(float* )mallocCUDADevice( sizeof(float) * blokline * tempDemCols);
d_demsloper_x=(float* )mallocCUDADevice( sizeof(float) * blokline * tempDemCols);
d_demsloper_y=(float* )mallocCUDADevice( sizeof(float) * blokline * tempDemCols);
d_demsloper_z=(float* )mallocCUDADevice( sizeof(float) * blokline * tempDemCols);
d_demsloper_angle= (float*)mallocCUDADevice( sizeof(float) * blokline * tempDemCols);
float* h_dem_theta; // <20><><EFBFBD>߷<EFBFBD><DFB7><EFBFBD>ͼ
float* h_dem_phi;
float* d_dem_theta;
float* d_dem_phi;
h_dem_theta=(float* )mallocCUDAHost( sizeof(float) * blokline * tempDemCols);
h_dem_phi= (float*)mallocCUDAHost(sizeof(float) * blokline * tempDemCols);
d_dem_theta= (float*)mallocCUDADevice( sizeof(float) * blokline * tempDemCols);// 9
d_dem_phi= (float*)mallocCUDADevice( sizeof(float) * blokline * tempDemCols);
HostToDevice((void*)h_dem_theta, (void*)d_dem_theta, sizeof(float) * blokline * tempDemCols);
HostToDevice((void*)h_dem_phi, (void*)d_dem_phi, sizeof(float) * blokline * tempDemCols);
// <20><>ǰ<EFBFBD><C7B0><EFBFBD><EFBFBD><EFBFBD><EFBFBD><EFBFBD><EFBFBD><EFBFBD><EFBFBD><EFBFBD><EFBFBD>
float* h_R;// <20><><EFBFBD><EFBFBD><E4B7BD>
float* h_localangle;//<2F><><EFBFBD><EFBFBD><EFBFBD><EFBFBD>
float* d_R;// <20><><EFBFBD><EFBFBD><E4B7BD>
float* d_localangle;//<2F><><EFBFBD><EFBFBD><EFBFBD><EFBFBD>
h_R=(float* )mallocCUDAHost(sizeof(float) * blokline * tempDemCols);
h_localangle= (float*)mallocCUDAHost(sizeof(float) * blokline * tempDemCols); // 11
d_R= (float*)mallocCUDADevice(sizeof(float) * blokline * tempDemCols);
d_localangle= (float*)mallocCUDADevice( sizeof(float) * blokline * tempDemCols);
float* h_RstX;
float* h_RstY;
float* h_RstZ;
float* d_RstX;
float* d_RstY;
float* d_RstZ;
h_RstX=(float*)mallocCUDAHost( sizeof(float) * blokline * tempDemCols);
h_RstY=(float*)mallocCUDAHost( sizeof(float) * blokline * tempDemCols);
h_RstZ=(float*)mallocCUDAHost( sizeof(float) * blokline * tempDemCols); // 14
d_RstX=(float*)mallocCUDADevice( sizeof(float) * blokline * tempDemCols);
d_RstY=(float*)mallocCUDADevice( sizeof(float) * blokline * tempDemCols);
d_RstZ=(float*)mallocCUDADevice( sizeof(float) * blokline * tempDemCols);
float* h_sigma0;
float* h_TransAnt;
float* h_ReciveAnt;
float* d_sigma0;
float* d_TransAnt;
float* d_ReciveAnt;
h_sigma0= (float*)mallocCUDAHost( sizeof(float)* blokline* tempDemCols);
h_TransAnt = (float*)mallocCUDAHost(sizeof(float) * blokline * tempDemCols);
h_ReciveAnt = (float*)mallocCUDAHost(sizeof(float) * blokline * tempDemCols); // 17
d_sigma0= (float*)mallocCUDADevice( sizeof(float) * blokline * tempDemCols);
d_TransAnt= (float*)mallocCUDADevice( sizeof(float) * blokline * tempDemCols);
d_ReciveAnt= (float*)mallocCUDADevice( sizeof(float) * blokline * tempDemCols);
// <20>ز<EFBFBD>
cuComplex* h_echoAmp;
cuComplex* d_echoAmp;
h_echoAmp=(cuComplex*)mallocCUDAHost(sizeof(cuComplex) * blokline * tempDemCols);
d_echoAmp=(cuComplex*)mallocCUDADevice( sizeof(cuComplex) * blokline * tempDemCols); //19
long* h_FreqID;
long* d_FreqID;
h_FreqID=(long*)mallocCUDAHost( sizeof(long) * blokline * tempDemCols);
d_FreqID=(long*)mallocCUDADevice( sizeof(long) * blokline * tempDemCols); //21
// <20>ر<EFBFBD><D8B1><EFBFBD><EFBFBD><EFBFBD><EFBFBD><EFBFBD><EFBFBD><EFBFBD>
Eigen::MatrixXd landcover = Eigen::MatrixXd::Zero(blokline, tempDemCols);// <20><><EFBFBD><EFBFBD><E6B8B2><EFBFBD><EFBFBD><EFBFBD><EFBFBD>
long* h_demcls = (long*)mallocCUDAHost(sizeof(long) * blokline * tempDemCols);
long* d_demcls = (long*)mallocCUDADevice(sizeof(long) * blokline * tempDemCols);
float* h_amp=(float*)mallocCUDAHost(sizeof(float) * blokline * tempDemCols);
float* d_amp=(float*)mallocCUDADevice(sizeof(float) * blokline * tempDemCols);
for (startline = 0; startline < demRow; startline = startline + blokline) {
long newblokline = blokline;
if ((startline + blokline) >= demRow) {
newblokline = demRow - startline;
bloklineflag = true;
}
dem_x = demxyz.getData(startline, 0, newblokline, demxyz.width, 1); // <20><><EFBFBD><EFBFBD><EFBFBD><EFBFBD><EFBFBD><EFBFBD>
dem_y = demxyz.getData(startline, 0, newblokline, demxyz.width, 2);
dem_z = demxyz.getData(startline, 0, newblokline, demxyz.width, 3);
demsloper_x = demsloperxyz.getData(startline, 0, newblokline, demsloperxyz.width, 1);
demsloper_y = demsloperxyz.getData(startline, 0, newblokline, demsloperxyz.width, 2);
demsloper_z = demsloperxyz.getData(startline, 0, newblokline, demsloperxyz.width, 3);
sloperAngle = demsloperxyz.getData(startline, 0, newblokline, demsloperxyz.width, 4);
landcover = demlandcls.getData(startline, 0, newblokline, demlandcls.width, 1);
if (bloklineflag) {
FreeCUDAHost(h_dem_x); FreeCUDADevice(d_dem_x);
FreeCUDAHost(h_dem_y); FreeCUDADevice(d_dem_y);
FreeCUDAHost(h_dem_z); FreeCUDADevice(d_dem_z);
FreeCUDAHost(h_demsloper_x); FreeCUDADevice(d_demsloper_x);
FreeCUDAHost(h_demsloper_y); FreeCUDADevice(d_demsloper_y);
FreeCUDAHost(h_demsloper_z); FreeCUDADevice(d_demsloper_z); //6
FreeCUDAHost(h_demsloper_angle); FreeCUDADevice(d_demsloper_angle);//7
FreeCUDAHost(h_dem_theta); FreeCUDADevice(d_dem_theta);
FreeCUDAHost(h_dem_phi); FreeCUDADevice(d_dem_phi); //9
FreeCUDAHost(h_R); FreeCUDADevice(d_R);
FreeCUDAHost(h_localangle); FreeCUDADevice(d_localangle); //11
FreeCUDAHost(h_RstX); FreeCUDADevice(d_RstX);
FreeCUDAHost(h_RstY); FreeCUDADevice(d_RstY);
FreeCUDAHost(h_RstZ); FreeCUDADevice(d_RstZ); //14
FreeCUDAHost(h_sigma0); FreeCUDADevice(d_sigma0);
FreeCUDAHost(h_TransAnt); FreeCUDADevice(d_TransAnt);
FreeCUDAHost(h_ReciveAnt); FreeCUDADevice(d_ReciveAnt); //17
FreeCUDAHost(h_echoAmp); FreeCUDADevice(d_echoAmp);//19
FreeCUDAHost(h_FreqID); FreeCUDADevice(d_FreqID);//20
FreeCUDAHost(h_demcls); FreeCUDADevice(d_demcls);
FreeCUDAHost(h_amp); FreeCUDADevice(d_amp);
h_dem_x = (float*)mallocCUDAHost(sizeof(float) * newblokline * tempDemCols);
h_dem_y = (float*)mallocCUDAHost(sizeof(float) * newblokline * tempDemCols);
h_dem_z = (float*)mallocCUDAHost(sizeof(float) * newblokline * tempDemCols);
h_demsloper_x = (float*)mallocCUDAHost(sizeof(float) * newblokline * tempDemCols);
h_demsloper_y = (float*)mallocCUDAHost(sizeof(float) * newblokline * tempDemCols);
h_demsloper_z = (float*)mallocCUDAHost(sizeof(float) * newblokline * tempDemCols);
h_demsloper_angle = (float*)mallocCUDAHost(sizeof(float) * blokline * tempDemCols);
h_dem_theta = (float*)mallocCUDAHost(sizeof(float) * newblokline * tempDemCols);
h_dem_phi = (float*)mallocCUDAHost(sizeof(float) * newblokline * tempDemCols);
h_R = (float*)mallocCUDAHost(sizeof(float) * newblokline * tempDemCols);
h_localangle = (float*)mallocCUDAHost(sizeof(float) * newblokline * tempDemCols);
h_RstX = (float*)mallocCUDAHost(sizeof(float) * newblokline * tempDemCols);
h_RstY = (float*)mallocCUDAHost(sizeof(float) * newblokline * tempDemCols);
h_RstZ = (float*)mallocCUDAHost(sizeof(float) * newblokline * tempDemCols);
h_sigma0 = (float*)mallocCUDAHost(sizeof(float) * newblokline * tempDemCols);
h_TransAnt = (float*)mallocCUDAHost(sizeof(float) * newblokline * tempDemCols);
h_ReciveAnt = (float*)mallocCUDAHost(sizeof(float) * newblokline * tempDemCols);
h_echoAmp = (cuComplex*)mallocCUDAHost(sizeof(cuComplex) * newblokline * tempDemCols);
h_FreqID = (long*)mallocCUDAHost(sizeof(long) * newblokline * tempDemCols);
h_demcls = (long*)mallocCUDAHost(sizeof(long) * newblokline * tempDemCols);
h_amp = (float*)mallocCUDAHost(sizeof(float) * newblokline * tempDemCols);
d_dem_x=(float*)mallocCUDADevice(sizeof(float) * newblokline * tempDemCols);
d_dem_y=(float*)mallocCUDADevice(sizeof(float) * newblokline * tempDemCols);
d_dem_z=(float*)mallocCUDADevice(sizeof(float) * newblokline * tempDemCols);
d_demsloper_x=(float*)mallocCUDADevice(sizeof(float) * newblokline * tempDemCols);
d_demsloper_y=(float*)mallocCUDADevice(sizeof(float) * newblokline * tempDemCols);
d_demsloper_z=(float*)mallocCUDADevice(sizeof(float) * newblokline * tempDemCols);//6
d_demsloper_angle=(float*)mallocCUDADevice(sizeof(float) * newblokline * tempDemCols);//7
d_dem_theta=(float*)mallocCUDADevice(sizeof(float) * newblokline * tempDemCols);
d_dem_phi=(float*)mallocCUDADevice(sizeof(float) * newblokline * tempDemCols);// 9
d_R=(float*)mallocCUDADevice(sizeof(float) * newblokline * tempDemCols);
d_localangle=(float*)mallocCUDADevice(sizeof(float) * newblokline * tempDemCols);
d_RstX=(float*)mallocCUDADevice(sizeof(float) * newblokline * tempDemCols);
d_RstY=(float*)mallocCUDADevice(sizeof(float) * newblokline * tempDemCols);
d_RstZ=(float*)mallocCUDADevice(sizeof(float) * newblokline * tempDemCols);
d_sigma0=(float*)mallocCUDADevice(sizeof(float) * newblokline * tempDemCols);
d_TransAnt=(float*)mallocCUDADevice(sizeof(float) * newblokline * tempDemCols);
d_ReciveAnt=(float*)mallocCUDADevice(sizeof(float) * newblokline * tempDemCols);
d_echoAmp=(cuComplex*)mallocCUDADevice(sizeof(cuComplex) * newblokline * tempDemCols);
d_FreqID=(long*)mallocCUDADevice(sizeof(long) * newblokline * tempDemCols);
d_demcls = (long*)mallocCUDADevice(sizeof(long) * newblokline * tempDemCols);
d_amp = (float*)mallocCUDADevice(sizeof(float) * newblokline * tempDemCols);
}
//# pragma omp parallel for
for (long i = 0; i < newblokline; i++) {
for (long j = 0; j < demxyz.width; j++) {
h_dem_x[i * demxyz.width + j] = float(dem_x(i, j));
h_dem_y[i * demxyz.width + j] = float(dem_y(i, j));
h_dem_z[i * demxyz.width + j] = float(dem_z(i, j));
h_demsloper_x[i * demxyz.width + j] = float(demsloper_x(i, j));
h_demsloper_y[i * demxyz.width + j] = float(demsloper_y(i, j));
h_demsloper_z[i * demxyz.width + j] = float(demsloper_z(i, j));
h_demsloper_angle[i * demxyz.width + j] = float(sloperAngle(i, j));
h_demcls[i * demxyz.width + j] = clamap[long(landcover(i, j))];
h_amp[i * demxyz.width + j] = 0.0f;
}
}
HostToDevice((void*)h_dem_x, (void*)d_dem_x, sizeof(float) * newblokline * tempDemCols); // <20><><EFBFBD><EFBFBD> <20><><EFBFBD><EFBFBD> -> GPU
HostToDevice((void*)h_dem_y, (void*)d_dem_y, sizeof(float) * newblokline * tempDemCols);
HostToDevice((void*)h_dem_z, (void*)d_dem_z, sizeof(float) * newblokline * tempDemCols);
HostToDevice((void*)h_demsloper_x, (void*)d_demsloper_x, sizeof(float) * newblokline * tempDemCols);
HostToDevice((void*)h_demsloper_y, (void*)d_demsloper_y, sizeof(float) * newblokline * tempDemCols);
HostToDevice((void*)h_demsloper_z, (void*)d_demsloper_z, sizeof(float) * newblokline * tempDemCols);
HostToDevice((void*)h_demsloper_angle, (void*)d_demsloper_angle, sizeof(float) * newblokline * tempDemCols);
HostToDevice((void*)h_dem_theta, (void*)d_dem_theta, sizeof(float) * newblokline * tempDemCols);
HostToDevice((void*)h_dem_phi, (void*)d_dem_phi, sizeof(float) * newblokline* tempDemCols);
HostToDevice((void*)h_demcls, (void*)d_demcls, sizeof(long) * newblokline* tempDemCols);
HostToDevice((void*)h_amp, (void*)d_amp, sizeof(float) * newblokline* tempDemCols);
#ifdef __PRFDEBUG__
testOutAmpArr("h_dem_x.bin", h_dem_x, newblokline, tempDemCols);
testOutAmpArr("h_dem_y.bin", h_dem_y, newblokline, tempDemCols);
testOutAmpArr("h_dem_z.bin", h_dem_z, newblokline, tempDemCols);
testOutClsArr( "h_demcls.bin" , h_demcls, newblokline, tempDemCols);
#endif // __PRFDEBUG__
long pixelcount = newblokline * tempDemCols;
long echoblockline = Memory1MB * 2000 / 8 / 2 / PlusePoint;
long startprfid = 0;
for (startprfid = 0; startprfid < pluseCount; startprfid = startprfid + echoblockline) {
long templine = startprfid + echoblockline < PluseCount ? echoblockline : PluseCount - startprfid;
std::shared_ptr<std::complex<double>> echotemp = this->EchoSimulationData->getEchoArr(startprfid, templine);
for (long tempprfid = 0; tempprfid < templine; tempprfid++) {
{// <20><><EFBFBD><EFBFBD>
long prfid = tempprfid + startprfid;
// <20><><EFBFBD><EFBFBD>λ<EFBFBD><CEBB>
float antpx = sateOirbtNodes[prfid].Px;
float antpy = sateOirbtNodes[prfid].Py;
float antpz = sateOirbtNodes[prfid].Pz;
float antvx = sateOirbtNodes[prfid].Vx;
float antvy = sateOirbtNodes[prfid].Vy;
float antvz = sateOirbtNodes[prfid].Vz; //6
float antdirectx = sateOirbtNodes[prfid].AntDirecX;
float antdirecty = sateOirbtNodes[prfid].AntDirecY;
float antdirectz = sateOirbtNodes[prfid].AntDirecZ; // 9 <20><><EFBFBD><EFBFBD>ָ<EFBFBD><D6B8>
float antXaxisX = sateOirbtNodes[prfid].AntXaxisX;
float antXaxisY = sateOirbtNodes[prfid].AntXaxisY;
float antXaxisZ = sateOirbtNodes[prfid].AntXaxisZ;//12 <20><><EFBFBD><EFBFBD><EFBFBD><EFBFBD><EFBFBD><EFBFBD>ϵ
float antYaxisX = sateOirbtNodes[prfid].AntYaxisX;
float antYaxisY = sateOirbtNodes[prfid].AntYaxisY;
float antYaxisZ = sateOirbtNodes[prfid].AntYaxisZ;//15
float antZaxisX = sateOirbtNodes[prfid].AntZaxisX;
float antZaxisY = sateOirbtNodes[prfid].AntZaxisY;
float antZaxisZ = sateOirbtNodes[prfid].AntZaxisZ;//18
#ifdef __PRFDEBUG__
std::cout << "ant Position=[" << antpx << "," << antpy << "," << antpz << "]" << std::endl;
#endif // __PRFDEBUG__
make_VectorA_B(antpx, antpy, antpz, d_dem_x, d_dem_y, d_dem_z, d_RstX, d_RstY, d_RstZ, pixelcount); // Rst = Rs - Rt; <20><><EFBFBD><EFBFBD>-> ָ<><D6B8>
Norm_Vector(d_RstX, d_RstY, d_RstZ, d_R, pixelcount); // R
cosAngle_VA_AB(d_RstX, d_RstY, d_RstZ, d_demsloper_x, d_demsloper_y, d_demsloper_z, d_localangle, pixelcount); // <20>ֲ<EFBFBD><D6B2><EFBFBD><EFBFBD><EFBFBD><EFBFBD><EFBFBD>
SatelliteAntDirectNormal(d_RstX, d_RstY, d_RstZ,
antXaxisX, antXaxisY, antXaxisZ,
antYaxisX, antYaxisY, antYaxisZ,
antZaxisX, antZaxisY, antZaxisZ,
antdirectx, antdirecty, antdirectz,
d_dem_theta, d_dem_phi, pixelcount);// <20><><EFBFBD><EFBFBD><EFBFBD>Ƕ<EFBFBD>
#ifdef __PRFDEBUG__
//DeviceToHost(h_RstX, d_RstX, sizeof(float)* pixelcount);
//DeviceToHost(h_RstY, d_RstY, sizeof(float)* pixelcount);
//DeviceToHost(h_RstZ, d_RstZ, sizeof(float)* pixelcount);
//testOutAmpArr("h_RstX.bin", h_RstX, newblokline, tempDemCols);
//testOutAmpArr("h_RstY.bin", h_RstY, newblokline, tempDemCols);
//testOutAmpArr("h_RstZ.bin", h_RstZ, newblokline, tempDemCols);
#endif // __PRFDEBUG__
AntPatternInterpGain(d_dem_theta, d_dem_phi, d_TransAnt, d_TantPattern, TstartTheta, TstartPhi, Tdtheta, Tdphi, Tthetanum, Tphinum, pixelcount);
AntPatternInterpGain(d_dem_theta, d_dem_phi, d_ReciveAnt, d_RantPattern, RstartTheta, RstartPhi, Rdtheta, Rdphi, Rthetanum, Rphinum, pixelcount);
#ifdef __PRFDEBUG__
DeviceToHost(h_dem_theta, d_dem_theta, sizeof(float)* pixelcount); // <20><>GPU -> <20><><EFBFBD><EFBFBD>
DeviceToHost(h_dem_phi, d_dem_phi, sizeof(float)* pixelcount);
DeviceToHost(h_localangle, d_localangle, sizeof(float)* pixelcount);
testOutAmpArr(QString("h_localangle_%1.bin").arg(prfid), h_localangle, newblokline, tempDemCols);
DeviceToHost(h_TransAnt, d_TransAnt, sizeof(float)* pixelcount);
DeviceToHost(h_ReciveAnt, d_ReciveAnt, sizeof(float)* pixelcount);
testOutAmpArr(QString("ant_theta_%1.bin").arg(prfid), h_dem_theta, newblokline, tempDemCols);
testOutAmpArr(QString("ant_phi_%1.bin").arg(prfid), h_dem_phi, newblokline, tempDemCols);
testOutAmpArr(QString("antPattern_Trans_%1.bin").arg(prfid), h_TransAnt, newblokline, tempDemCols);
testOutAmpArr(QString("antPattern_Receive_%1.bin").arg(prfid), h_ReciveAnt, newblokline, tempDemCols);
#endif // __PRFDEBUG__
CUDAInterpSigma(d_demcls, d_amp, d_localangle, pixelcount, d_clsSigmaParam, clamapid);
#ifdef __PRFDEBUG__
DeviceToHost(h_amp, d_amp, sizeof(float)* pixelcount);
testOutAmpArr(QString("amp_%1.bin").arg(prfid), h_amp,newblokline, tempDemCols);
#endif // __PRFDEBUG__
// <20><><EFBFBD><EFBFBD><EFBFBD>ز<EFBFBD>
calculationEcho(d_amp, d_TransAnt, d_ReciveAnt, d_localangle, d_R, d_demsloper_angle, NearRange, Fs, Pt, lamda, PlusePoint, d_echoAmp, d_FreqID, pixelcount);
DeviceToHost(h_echoAmp, d_echoAmp, sizeof(cuComplex) * pixelcount);
DeviceToHost(h_FreqID, d_FreqID, sizeof(long) * pixelcount);
//DeviceToHost(h_amp, d_amp, sizeof(float) * pixelcount);
#ifdef __PRFDEBUG__
float* h_echoAmp_real = (float*)mallocCUDAHost(sizeof(float) * pixelcount);
float* h_echoAmp_imag = (float*)mallocCUDAHost(sizeof(float) * pixelcount);
for (long freqi = 0; freqi < pixelcount; freqi++) {
h_echoAmp_real[freqi] = h_echoAmp[freqi].x;
h_echoAmp_imag[freqi] = h_echoAmp[freqi].y;
}
testOutAmpArr(QString("h_echoAmp_real_%1.bin").arg(prfid), h_echoAmp_real, newblokline, tempDemCols);
testOutAmpArr(QString("h_echoAmp_imag_%1.bin").arg(prfid), h_echoAmp_imag, newblokline, tempDemCols);
testOutClsArr(QString("h_FreqID_%1.bin").arg(prfid), h_FreqID, newblokline, tempDemCols);
FreeCUDAHost(h_echoAmp_real);
FreeCUDAHost(h_echoAmp_imag);
#endif // __PRFDEBUG__
for (long freqi = 0; freqi < pixelcount; freqi++) {
long pluseid = h_FreqID[freqi];
echotemp.get()[tempprfid * PlusePoint + pluseid] = std::complex<double>(h_echoAmp[freqi].x, h_echoAmp[freqi].y);
}
if (prfid % 1000 == 0) {
std::cout << "[" << QDateTime::currentDateTime().toString("yyyy-MM-dd hh:mm:ss.zzz").toStdString() << "] dem:\t" << startline << "\t-\t" << startline + newblokline << "\t:\t pluse :\t" << prfid << " / " << pluseCount << std::endl;
}
#ifdef __PRFDEBUG__
//this->EchoSimulationData->saveEchoArr(echotemp, startprfid, templine);
//exit(0);
#endif // __PRFDEBUG__
}
}
this->EchoSimulationData->saveEchoArr(echotemp, startprfid, templine);
}
}
std::cout << std::endl;
// <20><><EFBFBD><EFBFBD><EFBFBD><EFBFBD><EFBFBD><EFBFBD><EFBFBD>ͷ<EFBFBD>
FreeCUDAHost(h_dem_x); FreeCUDADevice(d_dem_x);
FreeCUDAHost(h_dem_y); FreeCUDADevice(d_dem_y);
FreeCUDAHost(h_dem_z); FreeCUDADevice(d_dem_z);
FreeCUDAHost(h_demsloper_x); FreeCUDADevice(d_demsloper_x);
FreeCUDAHost(h_demsloper_y); FreeCUDADevice(d_demsloper_y);
FreeCUDAHost(h_demsloper_z); FreeCUDADevice(d_demsloper_z); //6
FreeCUDAHost(h_demsloper_angle); FreeCUDADevice(d_demsloper_angle); //7
// <20><>ʱ<EFBFBD><CAB1><EFBFBD><EFBFBD><EFBFBD>ͷ<EFBFBD>
FreeCUDAHost(h_dem_theta); FreeCUDADevice(d_dem_theta);
FreeCUDAHost(h_dem_phi); FreeCUDADevice(d_dem_phi);// 9
FreeCUDAHost(h_R); FreeCUDADevice(d_R);
FreeCUDAHost(h_localangle); FreeCUDADevice(h_localangle); //11
FreeCUDAHost(h_RstX); FreeCUDADevice(d_RstX);
FreeCUDAHost(h_RstY); FreeCUDADevice(d_RstY);
FreeCUDAHost(h_RstZ); FreeCUDADevice(d_RstZ); //14
FreeCUDAHost(h_sigma0); FreeCUDADevice(d_sigma0);
FreeCUDAHost(h_TransAnt); FreeCUDADevice(d_TransAnt);
FreeCUDAHost(h_ReciveAnt); FreeCUDADevice(d_ReciveAnt); //17
FreeCUDAHost(h_echoAmp); FreeCUDADevice(d_echoAmp);//19
FreeCUDAHost(h_FreqID); FreeCUDADevice(d_FreqID);//20
FreeCUDAHost(h_demcls); FreeCUDADevice(d_demcls);
FreeCUDAHost(h_amp); FreeCUDADevice(d_amp);
#endif
this->EchoSimulationData->saveToXml();
return ErrorCode::SUCCESS;
}
void RTPCProcessMain(long num_thread, QString TansformPatternFilePath, QString ReceivePatternFilePath, QString simulationtaskName, QString OutEchoPath, QString GPSXmlPath, QString TaskXmlPath, QString demTiffPath, QString LandCoverPath, QString HHSigmaPath, QString HVSigmaPath, QString VHSigmaPath, QString VVSigmaPath)
{
std::shared_ptr < AbstractSARSatelliteModel> task = ReadSimulationSettingsXML(TaskXmlPath);
if (nullptr == task)
{
return;
}
else {
// <20><>ӡ<EFBFBD><D3A1><EFBFBD><EFBFBD>
qDebug() << "--------------Task Seting ---------------------------------------";
qDebug() << "SARImageStartTime: " << task->getSARImageStartTime();
qDebug() << "SARImageEndTime: " << task->getSARImageEndTime();
qDebug() << "BandWidth: " << task->getBandWidth();
qDebug() << "CenterFreq: " << task->getCenterFreq();
qDebug() << "PRF: " << task->getPRF();
qDebug() << "Fs: " << task->getFs();
qDebug() << "POLAR: " << task->getPolarType();
qDebug() << "NearRange: " << task->getNearRange();
qDebug() << "FarRange: " << task->getFarRange();
qDebug() << (task->getFarRange() - task->getNearRange()) * 2 / LIGHTSPEED * task->getFs();
qDebug() << "\n\n";
}
// 1.2 <20><><EFBFBD><EFBFBD><EFBFBD><EFBFBD><EFBFBD>߷<EFBFBD><DFB7><EFBFBD>ͼ
std::vector<RadiationPatternGainPoint> TansformPatternGainpoints = ReadGainFile(TansformPatternFilePath);
std::shared_ptr<AbstractRadiationPattern> TansformPatternGainPtr = CreateAbstractRadiationPattern(TansformPatternGainpoints);
std::vector<RadiationPatternGainPoint> ReceivePatternGainpoints = ReadGainFile(ReceivePatternFilePath);
std::shared_ptr<AbstractRadiationPattern> ReceivePatternGainPtr = CreateAbstractRadiationPattern(ReceivePatternGainpoints);
task->setTransformRadiationPattern(TansformPatternGainPtr);
task->setReceiveRadiationPattern(ReceivePatternGainPtr);
//2. <20><>ȡGPS<50>ڵ<EFBFBD>
std::vector<SatelliteOribtNode> nodes;
ErrorCode stateCode = ReadSateGPSPointsXML(GPSXmlPath, nodes);
if (stateCode != ErrorCode::SUCCESS)
{
qWarning() << QString::fromStdString(errorCode2errInfo(stateCode));
return;
}
else {}
std::shared_ptr<AbstractSatelliteOribtModel> SatelliteOribtModel = CreataPolyfitSatelliteOribtModel(nodes, task->getSARImageStartTime(), 3); // <20>Գ<EFBFBD><D4B3><EFBFBD><EFBFBD><EFBFBD>ʼʱ<CABC><CAB1><EFBFBD><EFBFBD>Ϊ ʱ<><CAB1><EFBFBD>ο<EFBFBD><CEBF><EFBFBD><EFBFBD><EFBFBD>
SatelliteOribtModel->setbeamAngle(task->getCenterLookAngle(), task->getIsRightLook()); // <20><><EFBFBD><EFBFBD><EFBFBD><EFBFBD><EFBFBD>߷<EFBFBD><DFB7><EFBFBD>ͼ
if (nullptr == SatelliteOribtModel)
{
return;
}
else {
task->setSatelliteOribtModel(SatelliteOribtModel);
}
qDebug() << "-------------- RTPC init ---------------------------------------";
RTPCProcessCls rtpc;
rtpc.setTaskSetting(task); //qDebug() << "setTaskSetting";
rtpc.setTaskFileName(simulationtaskName); //qDebug() << "setTaskFileName";
rtpc.setDEMTiffPath(demTiffPath); //qDebug() << "setDEMTiffPath";
rtpc.setLandCoverPath(LandCoverPath); //qDebug() << "setLandCoverPath";
rtpc.setHHSigmaPath(HHSigmaPath); //qDebug() << "setHHSigmaPath";
rtpc.setHVSigmaPath(HVSigmaPath); //qDebug() << "setHVSigmaPath";
rtpc.setVHSigmaPath(VHSigmaPath); //qDebug() << "setVHSigmaPath";
rtpc.setVVSigmaPath(VVSigmaPath); //qDebug() << "setVVSigmaPath";
rtpc.setOutEchoPath(OutEchoPath); //qDebug() << "setOutEchoPath";
qDebug() << "-------------- RTPC start---------------------------------------";
rtpc.Process(num_thread); // <20><><EFBFBD><EFBFBD><EFBFBD><EFBFBD><EFBFBD><EFBFBD>
qDebug() << "-------------- RTPC end---------------------------------------";
}
void testOutAntPatternTrans(QString antpatternfilename,float* antPatternArr,
double starttheta, double deltetheta,
double startphi, double deltaphi,
long thetanum, long phinum)
{
Eigen::MatrixXd antPatternMatrix(thetanum,phinum);
for (long t = 0; t < thetanum; ++t) {
for (long p = 0; p < phinum; ++p) {
long index = t * phinum + p;
if (index < thetanum * phinum) {
antPatternMatrix(t, p) = static_cast<double>(antPatternArr[index]); // Copy to Eigen matrix
}
}
}
Eigen::MatrixXd gt(2, 3);
gt(0, 0)=startphi;//x
gt(0, 1)=deltaphi;
gt(0, 2)=0;
gt(1, 0)=starttheta;
gt(1, 1)=0;
gt(1, 2)=deltetheta;
QString antpatternfilepath = getDebugDataPath(antpatternfilename);
gdalImage ds= CreategdalImage(antpatternfilepath, thetanum, phinum, 1, gt, "", true, true, true);
ds.saveImage(antPatternMatrix, 0, 0, 1);
}
void testOutClsArr(QString filename, long* amp, long rowcount, long colcount) {
Eigen::MatrixXd h_amp_img = Eigen::MatrixXd::Zero(rowcount, colcount);
for (long hii = 0; hii < rowcount; hii++) {
for (long hjj = 0; hjj < colcount; hjj++) {
h_amp_img(hii, hjj) = amp[hii * colcount + hjj];
}
}
QString ampPath = getDebugDataPath(filename);
saveEigenMatrixXd2Bin(h_amp_img, ampPath);
std::cout << filename.toLocal8Bit().constData() << std::endl;
std::cout << "max:\t" << h_amp_img.maxCoeff() << std::endl;
std::cout << "min:\t" << h_amp_img.minCoeff() << std::endl;
}
void testOutAmpArr(QString filename, float* amp, long rowcount, long colcount)
{
Eigen::MatrixXd h_amp_img = Eigen::MatrixXd::Zero(rowcount, colcount);
for (long hii = 0; hii < rowcount; hii++) {
for (long hjj = 0; hjj < colcount; hjj++) {
h_amp_img(hii, hjj) = amp[hii * colcount + hjj];
}
}
QString ampPath = getDebugDataPath(filename);
saveEigenMatrixXd2Bin(h_amp_img, ampPath);
std::cout << filename.toLocal8Bit().constData() << std::endl;
std::cout << "max:\t" << h_amp_img.maxCoeff() << std::endl;
std::cout << "min:\t" << h_amp_img.minCoeff() << std::endl;
}