RasterProcessTool/Toolbox/SimulationSARTool/SimulationSAR/RFPCProcessCls.cpp

1179 lines
46 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 <complex>
#include "GPUBaseTool.h"
#include "RFPCProcessCls.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"
#include "GPURFPC.cuh"
#include <cuda_runtime.h>
#include <cublas_v2.h>
#endif // __CUDANVCC___
//#include <Imageshow/ImageShowDialogClass.h>
CUDA_AntSate_PtrList* malloc_AntSate_PtrList(long PRFCount)
{
CUDA_AntSate_PtrList* antlist = (CUDA_AntSate_PtrList*)malloc(sizeof(CUDA_AntSate_PtrList));
antlist->h_antpx = (double*)mallocCUDAHost(sizeof(double) * PRFCount);
antlist->h_antpy = (double*)mallocCUDAHost(sizeof(double) * PRFCount);
antlist->h_antpz = (double*)mallocCUDAHost(sizeof(double) * PRFCount);
antlist->h_antvx = (double*)mallocCUDAHost(sizeof(double) * PRFCount);
antlist->h_antvy = (double*)mallocCUDAHost(sizeof(double) * PRFCount);
antlist->h_antvz = (double*)mallocCUDAHost(sizeof(double) * PRFCount);
antlist->h_antdirectx = (double*)mallocCUDAHost(sizeof(double) * PRFCount);
antlist->h_antdirecty = (double*)mallocCUDAHost(sizeof(double) * PRFCount);
antlist->h_antdirectz = (double*)mallocCUDAHost(sizeof(double) * PRFCount);
antlist->h_antXaxisX = (double*)mallocCUDAHost(sizeof(double) * PRFCount);
antlist->h_antXaxisY = (double*)mallocCUDAHost(sizeof(double) * PRFCount);
antlist->h_antXaxisZ = (double*)mallocCUDAHost(sizeof(double) * PRFCount);
antlist->h_antYaxisX = (double*)mallocCUDAHost(sizeof(double) * PRFCount);
antlist->h_antYaxisY = (double*)mallocCUDAHost(sizeof(double) * PRFCount);
antlist->h_antYaxisZ = (double*)mallocCUDAHost(sizeof(double) * PRFCount);
antlist->h_antZaxisX = (double*)mallocCUDAHost(sizeof(double) * PRFCount);
antlist->h_antZaxisY = (double*)mallocCUDAHost(sizeof(double) * PRFCount);
antlist->h_antZaxisZ = (double*)mallocCUDAHost(sizeof(double) * PRFCount);
antlist->d_antpx = (double*)mallocCUDADevice(sizeof(double) * PRFCount);
antlist->d_antpy = (double*)mallocCUDADevice(sizeof(double) * PRFCount);
antlist->d_antpz = (double*)mallocCUDADevice(sizeof(double) * PRFCount);
antlist->d_antvx = (double*)mallocCUDADevice(sizeof(double) * PRFCount);
antlist->d_antvy = (double*)mallocCUDADevice(sizeof(double) * PRFCount);
antlist->d_antvz = (double*)mallocCUDADevice(sizeof(double) * PRFCount);
antlist->d_antdirectx = (double*)mallocCUDADevice(sizeof(double) * PRFCount);
antlist->d_antdirecty = (double*)mallocCUDADevice(sizeof(double) * PRFCount);
antlist->d_antdirectz = (double*)mallocCUDADevice(sizeof(double) * PRFCount);
antlist->d_antXaxisX = (double*)mallocCUDADevice(sizeof(double) * PRFCount);
antlist->d_antXaxisY = (double*)mallocCUDADevice(sizeof(double) * PRFCount);
antlist->d_antXaxisZ = (double*)mallocCUDADevice(sizeof(double) * PRFCount);
antlist->d_antYaxisX = (double*)mallocCUDADevice(sizeof(double) * PRFCount);
antlist->d_antYaxisY = (double*)mallocCUDADevice(sizeof(double) * PRFCount);
antlist->d_antYaxisZ = (double*)mallocCUDADevice(sizeof(double) * PRFCount);
antlist->d_antZaxisX = (double*)mallocCUDADevice(sizeof(double) * PRFCount);
antlist->d_antZaxisY = (double*)mallocCUDADevice(sizeof(double) * PRFCount);
antlist->d_antZaxisZ = (double*)mallocCUDADevice(sizeof(double) * PRFCount);
antlist->PRF_len = PRFCount;
return antlist;
}
void Free_AntSate_PtrList(CUDA_AntSate_PtrList* antlist)
{
FreeCUDAHost(antlist->h_antpx);
FreeCUDAHost(antlist->h_antpy);
FreeCUDAHost(antlist->h_antpz);
FreeCUDAHost(antlist->h_antvx);
FreeCUDAHost(antlist->h_antvy);
FreeCUDAHost(antlist->h_antvz);
FreeCUDAHost(antlist->h_antdirectx);
FreeCUDAHost(antlist->h_antdirecty);
FreeCUDAHost(antlist->h_antdirectz);
FreeCUDAHost(antlist->h_antXaxisX);
FreeCUDAHost(antlist->h_antXaxisY);
FreeCUDAHost(antlist->h_antXaxisZ);
FreeCUDAHost(antlist->h_antYaxisX);
FreeCUDAHost(antlist->h_antYaxisY);
FreeCUDAHost(antlist->h_antYaxisZ);
FreeCUDAHost(antlist->h_antZaxisX);
FreeCUDAHost(antlist->h_antZaxisY);
FreeCUDAHost(antlist->h_antZaxisZ);
FreeCUDADevice(antlist->d_antpx);
FreeCUDADevice(antlist->d_antpy);
FreeCUDADevice(antlist->d_antpz);
FreeCUDADevice(antlist->d_antvx);
FreeCUDADevice(antlist->d_antvy);
FreeCUDADevice(antlist->d_antvz);
FreeCUDADevice(antlist->d_antdirectx);
FreeCUDADevice(antlist->d_antdirecty);
FreeCUDADevice(antlist->d_antdirectz);
FreeCUDADevice(antlist->d_antXaxisX);
FreeCUDADevice(antlist->d_antXaxisY);
FreeCUDADevice(antlist->d_antXaxisZ);
FreeCUDADevice(antlist->d_antYaxisX);
FreeCUDADevice(antlist->d_antYaxisY);
FreeCUDADevice(antlist->d_antYaxisZ);
FreeCUDADevice(antlist->d_antZaxisX);
FreeCUDADevice(antlist->d_antZaxisY);
FreeCUDADevice(antlist->d_antZaxisZ);
antlist->h_antpx = nullptr;
antlist->h_antpy = nullptr;
antlist->h_antpz = nullptr;
antlist->h_antvx = nullptr;
antlist->h_antvy = nullptr;
antlist->h_antvz = nullptr;
antlist->h_antdirectx = nullptr;
antlist->h_antdirecty = nullptr;
antlist->h_antdirectz = nullptr;
antlist->h_antXaxisX = nullptr;
antlist->h_antXaxisY = nullptr;
antlist->h_antXaxisZ = nullptr;
antlist->h_antYaxisX = nullptr;
antlist->h_antYaxisY = nullptr;
antlist->h_antYaxisZ = nullptr;
antlist->h_antZaxisX = nullptr;
antlist->h_antZaxisY = nullptr;
antlist->h_antZaxisZ = nullptr;
antlist->d_antpx = nullptr;
antlist->d_antpy = nullptr;
antlist->d_antpz = nullptr;
antlist->d_antvx = nullptr;
antlist->d_antvy = nullptr;
antlist->d_antvz = nullptr;
antlist->d_antdirectx = nullptr;
antlist->d_antdirecty = nullptr;
antlist->d_antdirectz = nullptr;
antlist->d_antXaxisX = nullptr;
antlist->d_antXaxisY = nullptr;
antlist->d_antXaxisZ = nullptr;
antlist->d_antYaxisX = nullptr;
antlist->d_antYaxisY = nullptr;
antlist->d_antYaxisZ = nullptr;
antlist->d_antZaxisX = nullptr;
antlist->d_antZaxisY = nullptr;
antlist->d_antZaxisZ = nullptr;
free(antlist);
antlist = nullptr;
}
void COPY_AntStation_FROM_HOST_GPU(std::shared_ptr<SatelliteOribtNode[]> sateOirbtNodes,
std::shared_ptr<CUDA_AntSate_PtrList> gpupptr,
long startPID,
long PRF_len)
{
assert(gpupptr->PRF_len <= PRF_len);
long prfid = 0;
for (long tempprfid = 0; tempprfid < PRF_len; tempprfid++) {
prfid = tempprfid + startPID;
gpupptr->h_antpx[tempprfid] = sateOirbtNodes[prfid].Px;
gpupptr->h_antpy[tempprfid] = sateOirbtNodes[prfid].Py;
gpupptr->h_antpz[tempprfid] = sateOirbtNodes[prfid].Pz;
gpupptr->h_antvx[tempprfid] = sateOirbtNodes[prfid].Vx;
gpupptr->h_antvy[tempprfid] = sateOirbtNodes[prfid].Vy;
gpupptr->h_antvz[tempprfid] = sateOirbtNodes[prfid].Vz; //6
gpupptr->h_antdirectx[tempprfid] = sateOirbtNodes[prfid].AntDirecX;
gpupptr->h_antdirecty[tempprfid] = sateOirbtNodes[prfid].AntDirecY;
gpupptr->h_antdirectz[tempprfid] = sateOirbtNodes[prfid].AntDirecZ;
gpupptr->h_antXaxisX[tempprfid] = sateOirbtNodes[prfid].AntXaxisX;
gpupptr->h_antXaxisY[tempprfid] = sateOirbtNodes[prfid].AntXaxisY;
gpupptr->h_antXaxisZ[tempprfid] = sateOirbtNodes[prfid].AntXaxisZ;//12
gpupptr->h_antYaxisX[tempprfid] = sateOirbtNodes[prfid].AntYaxisX;
gpupptr->h_antYaxisY[tempprfid] = sateOirbtNodes[prfid].AntYaxisY;
gpupptr->h_antYaxisZ[tempprfid] = sateOirbtNodes[prfid].AntYaxisZ;//15
gpupptr->h_antZaxisX[tempprfid] = sateOirbtNodes[prfid].AntZaxisX;
gpupptr->h_antZaxisY[tempprfid] = sateOirbtNodes[prfid].AntZaxisY;
gpupptr->h_antZaxisZ[tempprfid] = sateOirbtNodes[prfid].AntZaxisZ;//18
}
HostToDevice(gpupptr->h_antpx, gpupptr->d_antpx, sizeof(double) * PRF_len);
HostToDevice(gpupptr->h_antpy, gpupptr->d_antpy, sizeof(double) * PRF_len);
HostToDevice(gpupptr->h_antpz, gpupptr->d_antpz, sizeof(double) * PRF_len);
HostToDevice(gpupptr->h_antvx, gpupptr->d_antvx, sizeof(double) * PRF_len);
HostToDevice(gpupptr->h_antvy, gpupptr->d_antvy, sizeof(double) * PRF_len);
HostToDevice(gpupptr->h_antvz, gpupptr->d_antvz, sizeof(double) * PRF_len);
HostToDevice(gpupptr->h_antdirectx, gpupptr->d_antdirectx, sizeof(double) * PRF_len);
HostToDevice(gpupptr->h_antdirecty, gpupptr->d_antdirecty, sizeof(double) * PRF_len);
HostToDevice(gpupptr->h_antdirectz, gpupptr->d_antdirectz, sizeof(double) * PRF_len);
HostToDevice(gpupptr->h_antXaxisX, gpupptr->d_antXaxisX, sizeof(double) * PRF_len);
HostToDevice(gpupptr->h_antXaxisY, gpupptr->d_antXaxisY, sizeof(double) * PRF_len);
HostToDevice(gpupptr->h_antXaxisZ, gpupptr->d_antXaxisZ, sizeof(double) * PRF_len);
HostToDevice(gpupptr->h_antYaxisX, gpupptr->d_antYaxisX, sizeof(double) * PRF_len);
HostToDevice(gpupptr->h_antYaxisY, gpupptr->d_antYaxisY, sizeof(double) * PRF_len);
HostToDevice(gpupptr->h_antYaxisZ, gpupptr->d_antYaxisZ, sizeof(double) * PRF_len);
HostToDevice(gpupptr->h_antZaxisX, gpupptr->d_antZaxisX, sizeof(double) * PRF_len);
HostToDevice(gpupptr->h_antZaxisY, gpupptr->d_antZaxisY, sizeof(double) * PRF_len);
HostToDevice(gpupptr->h_antZaxisZ, gpupptr->d_antZaxisZ, sizeof(double) * PRF_len);
}
RFPCProcessCls::RFPCProcessCls()
{
this->PluseCount = 0;
this->PlusePoint = 0;
this->TaskSetting = nullptr;
this->EchoSimulationData = nullptr;
this->LandCoverPath = "";
this->OutEchoPath = "";
this->LandCoverPath.clear();
this->OutEchoPath.clear();
this->SigmaDatabasePtr = std::shared_ptr<SigmaDatabase>(new SigmaDatabase);
}
RFPCProcessCls::~RFPCProcessCls()
{
}
void RFPCProcessCls::setTaskSetting(std::shared_ptr < AbstractSARSatelliteModel> TaskSetting)
{
this->TaskSetting = std::shared_ptr < AbstractSARSatelliteModel>(TaskSetting);
qDebug() << "RFPCProcessCls::setTaskSetting";
}
void RFPCProcessCls::setEchoSimulationDataSetting(std::shared_ptr<EchoL0Dataset> EchoSimulationData)
{
this->EchoSimulationData = std::shared_ptr<EchoL0Dataset>(EchoSimulationData);
qDebug() << "RFPCProcessCls::setEchoSimulationDataSetting";
}
void RFPCProcessCls::setTaskFileName(QString EchoFileName)
{
this->TaskFileName = EchoFileName;
}
void RFPCProcessCls::setDEMTiffPath(QString DEMTiffPath)
{
this->demxyzPath = DEMTiffPath;
}
void RFPCProcessCls::setSloperPath(QString InSloperPath)
{
this->demsloperPath = InSloperPath;
}
void RFPCProcessCls::setLandCoverPath(QString LandCoverPath)
{
this->LandCoverPath = LandCoverPath;
}
void RFPCProcessCls::setOutEchoPath(QString OutEchoPath)
{
this->OutEchoPath = OutEchoPath;
}
ErrorCode RFPCProcessCls::Process(long num_thread)
{
// RFPC 算法
qDebug() << u8"params init ....";
ErrorCode stateCode = this->InitParams();
if (stateCode != ErrorCode::SUCCESS) {
return stateCode;
}
else {}
qDebug() << "RFPCMainProcess";
//return ErrorCode::SUCCESS;
stateCode = this->InitEchoMaskArray();
if (stateCode != ErrorCode::SUCCESS) {
return stateCode;
}
else {}
qDebug() << "InitEchoMaskArray";
//stateCode = this->RFPCMainProcess(num_thread);
// 初始化回波
this->EchoSimulationData->initEchoArr(std::complex<double>(0, 0));
//return ErrorCode::SUCCESS;
//stateCode = this->RFPCMainProcess_GPU();
stateCode = this->RFPCMainProcess_MultiGPU_NoAntPattern();
if (stateCode != ErrorCode::SUCCESS) {
return stateCode;
}
else {}
return ErrorCode::SUCCESS;
}
ErrorCode RFPCProcessCls::InitParams()
{
if (nullptr == this->TaskSetting || this->demxyzPath.isEmpty() ||
this->LandCoverPath.isEmpty() || this->demsloperPath.isEmpty()) {
return ErrorCode::RFPC_PARAMSISEMPTY;
}
else {
}
// 归一化绝对路径
this->OutEchoPath = QDir(this->OutEchoPath).absolutePath();
// 回波大小
double imgStart_end = this->TaskSetting->getSARImageEndTime() - this->TaskSetting->getSARImageStartTime();
this->PluseCount = ceil(imgStart_end * this->TaskSetting->getPRF());
// 远斜距
//double rangeTimeSample = * 2.0 / LIGHTSPEED;
double drange = LIGHTSPEED / 2.0 / this->TaskSetting->getBandWidth();
QVector<double> freqlist = this->TaskSetting->getFreqList();
const long freqnum = freqlist.count();
this->PlusePoint = freqnum;// ceil((this->TaskSetting->getFarRange() - this->TaskSetting->getNearRange()) / LIGHTSPEED * 2 * this->TaskSetting->getBandWidth());
this->TaskSetting->setFarRange(this->TaskSetting->getNearRange() + (this->PlusePoint-1) * drange);
//ceil(rangeTimeSample * this->TaskSetting->getFs());
// 初始化回波存放位置
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->setBandwidth(this->TaskSetting->getBandWidth());
this->EchoSimulationData->setCenterAngle(this->TaskSetting->getCenterLookAngle());
this->EchoSimulationData->setLookSide(this->TaskSetting->getIsRightLook() ? "R" : "L");
this->EchoSimulationData->OpenOrNew(OutEchoPath, TaskFileName, PluseCount, PlusePoint);
this->EchoSimulationData->setRefPhaseRange(this->TaskSetting->getRefphaseRange());
this->EchoSimulationData->saveToXml();
QString tmpfolderPath = QDir(OutEchoPath).filePath("tmp");
if (QDir(tmpfolderPath).exists() == false) {
QDir(OutEchoPath).mkpath(tmpfolderPath);
}
this->tmpfolderPath = tmpfolderPath;
return ErrorCode::SUCCESS;
}
ErrorCode RFPCProcessCls::InitEchoMaskArray()
{
QString name = this->EchoSimulationData->getSimulationTaskName();
this->OutEchoMaskPath = JoinPath(this->OutEchoPath, name + "_echomask.bin");
Eigen::MatrixXd gt(2, 3);
gt(0, 0) = 0;
gt(0, 1) = 1;
gt(0, 2) = 0;
gt(1, 0) = 0;
gt(1, 1) = 0;
gt(1, 2) = 1;
gdalImage echomaskImg = CreategdalImage(this->OutEchoMaskPath,
this->EchoSimulationData->getPluseCount(),
this->EchoSimulationData->getPlusePoints(),
1,
gt, "",
false, true, true);
long cols = this->EchoSimulationData->getPlusePoints();
long rows = this->EchoSimulationData->getPluseCount();
long blocksize = Memory1GB / 8 / this->EchoSimulationData->getPlusePoints() * 4;
for (long startid = 0; startid < rows; startid = startid + blocksize) {
Eigen::MatrixXd data = echomaskImg.getData(startid, 0, blocksize, cols, 1);
data = data.array() * 0;
echomaskImg.saveImage(data, startid, 0, 1);
}
return ErrorCode::SUCCESS;
}
std::shared_ptr<SatelliteOribtNode[]> RFPCProcessCls::getSatelliteOribtNodes(double prf_time, double dt, bool antflag, long double imageStarttime)
{
std::shared_ptr<SatelliteOribtNode[]> sateOirbtNodes(new SatelliteOribtNode[this->PluseCount], delArrPtr);
{ // 姿态计算不同
qDebug() << "Ant position finished started !!!";
// 计算姿态
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; // 加速度
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 !!!";
}
return sateOirbtNodes;
}
void RFPCProcessMain(long num_thread,
QString TansformPatternFilePath, QString ReceivePatternFilePath,
QString simulationtaskName, QString OutEchoPath,
QString GPSXmlPath, QString TaskXmlPath,QString demTiffPath, QString sloperPath, QString LandCoverPath)
{
std::shared_ptr < AbstractSARSatelliteModel> task = ReadSimulationSettingsXML(TaskXmlPath);
if (nullptr == task)
{
return;
}
else {
// 打印参数
qDebug() << "--------------Task Seting ---------------------------------------";
qDebug() << "SARImageStartTime: " << task->getSARImageStartTime();
qDebug() << "SARImageEndTime: " << task->getSARImageEndTime();
qDebug() << "BandWidth: " << task->getBandWidth();
qDebug() << "CenterFreq: " << task->getCenterFreq();
qDebug() << "PRF: " << task->getPRF();
qDebug() << "Lamda: " << task->getCenterLamda();
qDebug() << "Fs: " << task->getFs();
qDebug() << "POLAR: " << task->getPolarType();
qDebug() << "NearRange: " << task->getNearRange();
qDebug() << "FarRange: " << task->getFarRange();
qDebug() << "RefRange: " << task->getRefphaseRange();
qDebug() << (task->getFarRange() - task->getNearRange()) * 2 / LIGHTSPEED * task->getFs();
qDebug() << "\n\n";
}
// 1.2 设置天线方向图
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. 读取GPS节点
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); // 以成像开始时间作为 时间参考起点
SatelliteOribtModel->setbeamAngle(task->getCenterLookAngle(), task->getIsRightLook()); // 设置天线方向图
if (nullptr == SatelliteOribtModel)
{
return;
}
else {
task->setSatelliteOribtModel(SatelliteOribtModel);
}
qDebug() << "-------------- RFPC init ---------------------------------------";
RFPCProcessCls RFPC;
RFPC.setTaskSetting(task); //qDebug() << "setTaskSetting";
RFPC.setTaskFileName(simulationtaskName); //qDebug() << "setTaskFileName";
RFPC.setDEMTiffPath(demTiffPath); //qDebug() << "setDEMTiffPath";
RFPC.setSloperPath(sloperPath); //qDebug() << "setDEMTiffPath";
RFPC.setLandCoverPath(LandCoverPath); //qDebug() << "setLandCoverPath";
RFPC.setOutEchoPath(OutEchoPath); //qDebug() << "setOutEchoPath";
qDebug() << "-------------- RFPC start---------------------------------------";
RFPC.Process(num_thread); // 处理程序
qDebug() << "-------------- RFPC end---------------------------------------";
}
ErrorCode RFPCProcessCls::RFPCMainProcess_GPU() {
/** 内存分配***************************************************/
long TargetMemoryMB = 500;
/** 参数区域***************************************************/
QVector<double> freqlist = this->TaskSetting->getFreqList();
const long freqnum = freqlist.count();
long freqnum_temp = freqnum;
float f0 = float(freqlist[0] / 1e9);
float dfreq = float((freqlist[1] - freqlist[0]) / 1e9);
long PRFCount = this->EchoSimulationData->getPluseCount();
double NearRange = this->EchoSimulationData->getNearRange(); // 近斜距
double FarRange = this->EchoSimulationData->getFarRange();
double Pt = this->TaskSetting->getPt() * this->TaskSetting->getGri();// 发射电压 1v
double lamda = this->TaskSetting->getCenterLamda(); // 波长
double refphaseRange = this->TaskSetting->getRefphaseRange(); // 参考相位斜距
double prf_time = 0;
double dt = 1 / this->TaskSetting->getPRF();// 获取每次脉冲的时间间隔
bool antflag = true; // 计算天线方向图
long double imageStarttime = this->TaskSetting->getSARImageStartTime();
qDebug() << "freqnum: " << freqnum << " f0: " << f0 << " dfreq: " << dfreq << "freqnum_temp: " << freqnum_temp;
this->EchoSimulationData->getAntPos();
std::shared_ptr<SatelliteOribtNode[]> sateOirbtNodes = this->getSatelliteOribtNodes(prf_time, dt, antflag, imageStarttime);
/** 天线方向图***************************************************/
std::shared_ptr<AbstractRadiationPattern> TransformPattern = this->TaskSetting->getTransformRadiationPattern(); // 发射天线方向图
std::shared_ptr<AbstractRadiationPattern> ReceivePattern = this->TaskSetting->getReceiveRadiationPattern(); // 接收天线方向图
POLARTYPEENUM polartype = this->TaskSetting->getPolarType();
PatternImageDesc TantPatternDesc = {};
double* h_TantPattern = nullptr;
double* d_TantPattern = nullptr;
double maxTransAntPatternValue = 0;
{
// 处理发射天线方向图
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);
h_TantPattern = (double*)mallocCUDAHost(sizeof(double) * Tthetanum * Tphinum);
d_TantPattern = (double*)mallocCUDADevice(sizeof(double) * 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);
maxTransAntPatternValue = powf(10.0, h_TantPattern[0] / 10);
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); // 转换为线性值
if (maxTransAntPatternValue < h_TantPattern[i * Tphinum + j]) {
maxTransAntPatternValue = h_TantPattern[i * Tphinum + j];
}
}
}
HostToDevice(h_TantPattern, d_TantPattern, sizeof(double) * Tthetanum * Tphinum);
TantPatternDesc.startTheta = TstartTheta;
TantPatternDesc.startPhi = TstartPhi;
TantPatternDesc.dtheta = Tdtheta;
TantPatternDesc.dphi = Tdphi;
TantPatternDesc.phinum = Tphinum;
TantPatternDesc.thetanum = Tthetanum;
}
PatternImageDesc RantPatternDesc = {};
double* h_RantPattern = nullptr;
double* d_RantPattern = nullptr;
double maxReceiveAntPatternValue = 0;
{
// 处理接收天线方向图
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);
h_RantPattern = (double*)mallocCUDAHost(sizeof(double) * Rthetanum * Rphinum);
d_RantPattern = (double*)mallocCUDADevice(sizeof(double) * 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);
maxReceiveAntPatternValue = powf(10.0, h_RantPattern[0] / 10);
for (long i = 0; i < Rthetanum; i++) {
for (long j = 0; j < Rphinum; j++) {
h_RantPattern[i * Rphinum + j] = powf(10.0, h_RantPattern[i * Rphinum + j] / 10);
if (maxReceiveAntPatternValue < h_RantPattern[i * Rphinum + j]) {
maxReceiveAntPatternValue = h_RantPattern[i * Rphinum + j];
}
}
}
HostToDevice(h_RantPattern, d_RantPattern, sizeof(double) * Rthetanum * Rphinum);
RantPatternDesc.startTheta = RstartTheta;
RantPatternDesc.startPhi = RstartPhi;
RantPatternDesc.dtheta = Rdtheta;
RantPatternDesc.dphi = Rdphi;
RantPatternDesc.phinum = Rphinum;
RantPatternDesc.thetanum = Rthetanum;
}
/** 坐标区域点***************************************************/
gdalImage demxyz(this->demxyzPath);// 地面点坐标
gdalImage demlandcls(this->LandCoverPath);// 地表覆盖类型
gdalImage demsloperxyz(this->demsloperPath);// 地面坡向
long demRow = demxyz.height;
long demCol = demxyz.width;
//处理地表覆盖
QMap<long, long> clamap;
long clamapid = 0;
long startline = 0;
{
long blokline = getBlockRows(2e4, demCol, sizeof(double),demRow);
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;
}
}
}
}
qDebug() << "class id recoding" ;
for (long id : clamap.keys()) {
qDebug() << id << " -> " << clamap[id] ;
}
}
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;
}
// 打印日志
qDebug() << "sigma params:" ;
qDebug() << "classid:\tp1\tp2\tp3\tp4\tp5\tp6" ;
for (long ii = 0; ii < clamapid; ii++) {
qDebug() << ii << ":\t" << h_clsSigmaParam[ii].p1;
qDebug() << "\t" << h_clsSigmaParam[ii].p2;
qDebug() << "\t" << h_clsSigmaParam[ii].p3;
qDebug() << "\t" << h_clsSigmaParam[ii].p4;
qDebug() << "\t" << h_clsSigmaParam[ii].p5;
qDebug() << "\t" << h_clsSigmaParam[ii].p6 ;
}
qDebug() << "";
}
HostToDevice(h_clsSigmaParam, d_clsSigmaParam, sizeof(CUDASigmaParam) * clamapid);
qDebug() << "CUDA class Proces finished!!!";
// 处理地面坐标
long blockline = getBlockRows(TargetMemoryMB, demCol, sizeof(double), demRow);
double* h_dem_x = (double*)mallocCUDAHost(sizeof(double) * blockline * demCol);
double* h_dem_y = (double*)mallocCUDAHost(sizeof(double) * blockline * demCol);
double* h_dem_z = (double*)mallocCUDAHost(sizeof(double) * blockline * demCol);
double* h_demsloper_x = (double*)mallocCUDAHost(sizeof(double) * blockline * demCol);
double* h_demsloper_y = (double*)mallocCUDAHost(sizeof(double) * blockline * demCol);
double* h_demsloper_z = (double*)mallocCUDAHost(sizeof(double) * blockline * demCol);
long* h_demcls = (long*)mallocCUDAHost(sizeof(long) * blockline * demCol);
double* d_dem_x = (double*)mallocCUDADevice(sizeof(double) * blockline * demCol);
double* d_dem_y = (double*)mallocCUDADevice(sizeof(double) * blockline * demCol);
double* d_dem_z = (double*)mallocCUDADevice(sizeof(double) * blockline * demCol);
double* d_demsloper_x = (double*)mallocCUDADevice(sizeof(double) * blockline * demCol);
double* d_demsloper_y = (double*)mallocCUDADevice(sizeof(double) * blockline * demCol);
double* d_demsloper_z = (double*)mallocCUDADevice(sizeof(double) * blockline * demCol);
long* d_demcls = (long*) mallocCUDADevice(sizeof(long) * blockline * demCol);
/** 处理回波***************************************************/
long echo_block_rows = getBlockRows(1000, freqnum, sizeof(float)*2, PRFCount);
float* h_echo_block_real = (float*)mallocCUDAHost(sizeof(float) * echo_block_rows * freqnum);
float* h_echo_block_imag = (float*)mallocCUDAHost(sizeof(float) * echo_block_rows * freqnum);
float* d_echo_block_real = (float*)mallocCUDADevice(sizeof(float) * echo_block_rows * freqnum);
float* d_echo_block_imag = (float*)mallocCUDADevice(sizeof(float) * echo_block_rows * freqnum);
float* d_temp_R = (float*)mallocCUDADevice(sizeof(float) * echo_block_rows * SHAREMEMORY_FLOAT_HALF); //2GB 距离
float* d_temp_amp = (float*)mallocCUDADevice(sizeof(float) * echo_block_rows * SHAREMEMORY_FLOAT_HALF);//2GB 强度
/** 主流程处理 ***************************************************/
qDebug() << "CUDA Main Proces";
for (long sprfid = 0; sprfid < PRFCount; sprfid = sprfid + echo_block_rows) {
long PRF_len = (sprfid + echo_block_rows) < PRFCount ? echo_block_rows : (PRFCount - sprfid);
qDebug() << "Start PRF: " << sprfid << "\t-\t" << sprfid + PRF_len << "\t:copy ant list host -> GPU";
std::shared_ptr< CUDA_AntSate_PtrList> antptrlist(malloc_AntSate_PtrList(PRF_len), Free_AntSate_PtrList);
COPY_AntStation_FROM_HOST_GPU(sateOirbtNodes, antptrlist, sprfid, PRF_len);
qDebug() << "freqnum: " << freqnum << " f0: " << f0 << " dfreq: " << dfreq << "freqnum_temp: " << freqnum_temp;
qDebug() << "Start PRF: " << sprfid << "\t-\t" << sprfid + PRF_len << "\t:copy echo data list host -> GPU";
std::shared_ptr<std::complex<double>> echo_temp = this->EchoSimulationData->getEchoArr(sprfid, PRF_len);
for (long ii = 0; ii < PRF_len; ii++) {
for (long jj = 0; jj < freqnum; jj++) {
h_echo_block_real[ii * freqnum + jj]=echo_temp.get()[ii * freqnum + jj].real();
h_echo_block_imag[ii * freqnum + jj]=echo_temp.get()[ii * freqnum + jj].imag();
}
}
HostToDevice(h_echo_block_real, d_echo_block_real, sizeof(float) * PRF_len* freqnum);
HostToDevice(h_echo_block_imag, d_echo_block_imag, sizeof(float) * PRF_len* freqnum);
for (startline = 0; startline < demRow; startline = startline + blockline) {
Eigen::MatrixXd dem_x = demxyz.getData(startline, 0, blockline, demCol, 1); // 地面坐标
Eigen::MatrixXd dem_y = demxyz.getData(startline, 0, blockline, demCol, 2);
Eigen::MatrixXd dem_z = demxyz.getData(startline, 0, blockline, demCol, 3);
Eigen::MatrixXd demsloper_x = demsloperxyz.getData(startline, 0, blockline, demCol, 1);
Eigen::MatrixXd demsloper_y = demsloperxyz.getData(startline, 0, blockline, demCol, 2);
Eigen::MatrixXd demsloper_z = demsloperxyz.getData(startline, 0, blockline, demCol, 3);
Eigen::MatrixXd landcover = demlandcls.getData(startline, 0, blockline, demCol, 1);
long temp_dem_row = dem_x.rows();
long temp_dem_col = dem_x.cols();
long temp_dem_count = dem_x.count();
// 更新数据格式
for (long i = 0; i < temp_dem_row; i++) {
for (long j = 0; j < temp_dem_col; j++) {
h_dem_x[i * temp_dem_col + j] = double(dem_x(i, j));
h_dem_y[i * temp_dem_col + j] = double(dem_y(i, j));
h_dem_z[i * temp_dem_col + j] = double(dem_z(i, j));
h_demsloper_x[i * temp_dem_col + j] = double(demsloper_x(i, j));
h_demsloper_y[i * temp_dem_col + j] = double(demsloper_y(i, j));
h_demsloper_z[i * temp_dem_col + j] = double(demsloper_z(i, j));
h_demcls[i * temp_dem_col + j] = clamap[long(landcover(i, j))];
}
}
qDebug() << "Start PRF: " << sprfid << "\t-\t" << sprfid + PRF_len << "\t:copy target data ("<< startline<<" - "<< startline + blockline << ") host -> GPU";
HostToDevice(h_dem_x, d_dem_x , sizeof(double) * blockline * demCol);
HostToDevice(h_dem_y, d_dem_y , sizeof(double) * blockline * demCol);
HostToDevice(h_dem_z, d_dem_z , sizeof(double) * blockline * demCol);
HostToDevice(h_demsloper_x, d_demsloper_x , sizeof(double) * blockline * demCol);
HostToDevice(h_demsloper_y, d_demsloper_y , sizeof(double) * blockline * demCol);
HostToDevice(h_demsloper_z, d_demsloper_z , sizeof(double) * blockline * demCol);
HostToDevice(h_demcls, d_demcls ,sizeof(long)* blockline* demCol);
// 分块处理
qDebug() << "Start PRF: " << sprfid << "\t-\t" << sprfid + PRF_len << "\t:GPU Computer target data (" << startline << "-" << startline + blockline << ")";
CUDA_RFPC_MainProcess(
antptrlist->d_antpx, antptrlist->d_antpy, antptrlist->d_antpz,
antptrlist->d_antXaxisX, antptrlist->d_antXaxisY, antptrlist->d_antXaxisZ, // 天线坐标系的X轴
antptrlist->d_antYaxisX, antptrlist->d_antYaxisY, antptrlist->d_antYaxisZ,// 天线坐标系的Y轴
antptrlist->d_antZaxisX, antptrlist->d_antZaxisY, antptrlist->d_antZaxisZ,// 天线坐标系的Z轴
antptrlist->d_antdirectx, antptrlist->d_antdirecty, antptrlist->d_antdirectz,// 天线的指向
PRF_len, freqnum,
f0,dfreq,
Pt,
refphaseRange,
// 天线方向图
d_TantPattern,
TantPatternDesc.startTheta, TantPatternDesc.startPhi, TantPatternDesc.dtheta, TantPatternDesc.dphi, TantPatternDesc.thetanum, TantPatternDesc.phinum,
d_RantPattern,
RantPatternDesc.startTheta, RantPatternDesc.startPhi, RantPatternDesc.dtheta, RantPatternDesc.dphi, RantPatternDesc.thetanum, RantPatternDesc.phinum,
maxTransAntPatternValue, maxReceiveAntPatternValue,
NearRange, FarRange, // 近斜据
d_dem_x, d_dem_y, d_dem_z, d_demcls, temp_dem_count, // 地面坐标
d_demsloper_x, d_demsloper_y, d_demsloper_z, // 地表坡度矢量
d_clsSigmaParam, clamapid,
d_echo_block_real, d_echo_block_imag,// 输出回波
d_temp_R, d_temp_amp
);
PRINT("dem : %d ~ %d / %d , echo: %d ~ %d / %d \n", startline, startline+ temp_dem_row, demRow, sprfid, sprfid+ PRF_len, PRFCount);
}
#if (defined __PRFDEBUG__) && (defined __PRFDEBUG_PRFINF__)
float* h_temp_R = (float*)mallocCUDAHost(sizeof(float) * echo_block_rows * SHAREMEMORY_FLOAT_HALF); //2GB 距离
float* h_temp_amp = (float*)mallocCUDAHost(sizeof(float) * echo_block_rows * SHAREMEMORY_FLOAT_HALF);//2GB 强度
DeviceToHost(h_temp_R, d_temp_R, sizeof(float) * echo_block_rows * SHAREMEMORY_FLOAT_HALF);
DeviceToHost(h_temp_amp, d_temp_amp, sizeof(float) * echo_block_rows * SHAREMEMORY_FLOAT_HALF);
testOutAmpArr("temp_R.bin", h_temp_R, echo_block_rows, SHAREMEMORY_FLOAT_HALF);
testOutAmpArr("temp_Amp.bin", h_temp_amp, echo_block_rows, SHAREMEMORY_FLOAT_HALF);
FreeCUDAHost(h_temp_R);
FreeCUDAHost(h_temp_amp);
#endif
DeviceToHost(h_echo_block_real, d_echo_block_real, sizeof(float) * PRF_len * freqnum);
DeviceToHost(h_echo_block_imag, d_echo_block_imag, sizeof(float) * PRF_len * freqnum);
qDebug() << "freqnum: " << freqnum << " f0: " << f0 << " dfreq: " << dfreq;
for (long ii = 0; ii < PRF_len; ii++) {
for (long jj = 0; jj < freqnum; jj++) {
echo_temp.get()[ii * freqnum + jj].real(h_echo_block_real[ii * freqnum + jj]);
echo_temp.get()[ii * freqnum + jj].imag(h_echo_block_imag[ii * freqnum + jj]);
}
}
this->EchoSimulationData->saveEchoArr(echo_temp, sprfid, PRF_len);
}
/** 内存释放***************************************************/
FreeCUDAHost(h_TantPattern);
FreeCUDAHost(h_RantPattern);
FreeCUDADevice(d_TantPattern);
FreeCUDADevice(d_RantPattern);
FreeCUDAHost(h_dem_x);
FreeCUDAHost(h_dem_y);
FreeCUDAHost(h_dem_z);
FreeCUDAHost(h_demsloper_x);
FreeCUDAHost(h_demsloper_y);
FreeCUDAHost(h_demsloper_z);
FreeCUDAHost(h_demcls);
FreeCUDAHost(h_echo_block_real);
FreeCUDAHost(h_echo_block_imag);
FreeCUDADevice(d_dem_x);
FreeCUDADevice(d_dem_y);
FreeCUDADevice(d_dem_z);
FreeCUDADevice(d_demsloper_x);
FreeCUDADevice(d_demsloper_y);
FreeCUDADevice(d_demsloper_z);
FreeCUDADevice(d_demcls);
FreeCUDADevice(d_echo_block_real);
FreeCUDADevice(d_echo_block_imag);
FreeCUDADevice(d_temp_R);
FreeCUDADevice(d_temp_amp);
return ErrorCode::SUCCESS;
}
ErrorCode RFPCProcessCls::RFPCMainProcess_MultiGPU_NoAntPattern()
{
int num_devices=0;
cudaGetDeviceCount(&num_devices);
PRINT("GPU Count : %d \n", num_devices);
long prfcount = this->EchoSimulationData->getPluseCount();
size_t prfblockcount = (prfcount + num_devices +2- 1) / num_devices;
PRINT("PRF COUNT : %d , child PRF COUNT: %d\n", prfcount, prfblockcount);
double prf_time = 0;
double dt = 1 / this->TaskSetting->getPRF();// 获取每次脉冲的时间间隔
bool antflag = true; // 计算天线方向图
long double imageStarttime = this->TaskSetting->getSARImageStartTime();
std::shared_ptr<SatelliteOribtNode[]> sateOirbtNodes = this->getSatelliteOribtNodes(prf_time, dt, antflag, imageStarttime);
#pragma omp parallel for
for (int devid = 0; devid < num_devices; devid++) {
cudaSetDevice(devid); // 确保当前线程操作指定的GPU设备
size_t startTid = devid * prfblockcount;
size_t prf_devLen = prfblockcount;
prf_devLen = (startTid + prf_devLen) < prfcount ? prf_devLen : (prfcount - startTid);
PRINT("dev ID:%d,start PRF ID: %d , PRF COUNT: %d \n", devid, startTid, prf_devLen);
this->RFPCMainProcess_GPU_NoAntPattern(startTid, prf_devLen, devid);
}
return ErrorCode::SUCCESS;
}
ErrorCode RFPCProcessCls::RFPCMainProcess_GPU_NoAntPattern(size_t startprfid, size_t prfcount, int devId)
{
PRINT("dev ID:%d,start PRF ID: %d , PRF COUNT: %d \n", devId,startprfid,prfcount);
/// 显存不限制
cudaSetDevice(devId); // 确保当前线程操作指定的GPU设备
POLARTYPEENUM polartype = this->TaskSetting->getPolarType();
std::map<long, SigmaParam> clssigmaParamsDict = this->SigmaDatabasePtr->getsigmaParams(polartype);;
std::map<long, CUDASigmaParam> clsCUDASigmaParamsDict;
for (const auto& pair : clssigmaParamsDict) {
clsCUDASigmaParamsDict.insert(std::pair<long, CUDASigmaParam>(pair.first,
CUDASigmaParam{
pair.second.p1,
pair.second.p2,
pair.second.p3,
pair.second.p4,
pair.second.p5,
pair.second.p6
}));
printf("clsid:%d, params: %e,%e,%e,%e,%e,%e \n", pair.first,
pair.second.p1,
pair.second.p2,
pair.second.p3,
pair.second.p4,
pair.second.p5,
pair.second.p6
);
}
// 读取类别
gdalImage demxyz(this->demxyzPath);// 地面点坐标
gdalImage demlandcls(this->LandCoverPath);// 地表覆盖类型
gdalImage slpxyz(this->demsloperPath);// 地面坡向
// 处理地面坐标
long demRow = demxyz.height;
long demCol = demxyz.width;
size_t demCount = size_t(demRow) * size_t(demCol);
std::shared_ptr<double> demX = readDataArr<double>(demxyz, 0, 0, demRow, demCol, 1, GDALREADARRCOPYMETHOD::VARIABLEMETHOD);
std::shared_ptr<double> demY = readDataArr<double>(demxyz, 0, 0, demRow, demCol, 2, GDALREADARRCOPYMETHOD::VARIABLEMETHOD);
std::shared_ptr<double> demZ = readDataArr<double>(demxyz, 0, 0, demRow, demCol, 3, GDALREADARRCOPYMETHOD::VARIABLEMETHOD);
std::shared_ptr<double> slpX = readDataArr<double>(slpxyz, 0, 0, demRow, demCol, 1, GDALREADARRCOPYMETHOD::VARIABLEMETHOD);
std::shared_ptr<double> slpY = readDataArr<double>(slpxyz, 0, 0, demRow, demCol, 2, GDALREADARRCOPYMETHOD::VARIABLEMETHOD);
std::shared_ptr<double> slpZ = readDataArr<double>(slpxyz, 0, 0, demRow, demCol, 3, GDALREADARRCOPYMETHOD::VARIABLEMETHOD);
std::shared_ptr<long> clsArr = readDataArr<long>(demlandcls, 0, 0, demRow, demCol, 1, GDALREADARRCOPYMETHOD::VARIABLEMETHOD);
// 检索类别数量
std::map<long, size_t> clsCountDict;
for (const auto& pair : clssigmaParamsDict) {
clsCountDict.insert(std::pair<long, size_t>(pair.first, 0));
}
for (size_t i = 0; i < demCount; i++) {
long clsid = clsArr.get()[i];
if (clsCountDict.find(clsid) != clsCountDict.end()) {
clsCountDict[clsid] = clsCountDict[clsid] + 1;
}
}
std::map<long, std::shared_ptr<GoalState>> clsGoalStateDict;
for (const auto& pair : clsCountDict) {
if (pair.second > 0) {
clsGoalStateDict.insert(
std::pair<long, std::shared_ptr<GoalState>>(
pair.first,
std::shared_ptr<GoalState>((GoalState*)mallocCUDAHost(sizeof(GoalState) * pair.second), FreeCUDAHost)));
PRINT("clsid : %d ,Count: %d\n", pair.first, pair.second);
}
}
// 分块处理大小
size_t blocksize = 1000;
std::map<long, size_t> clsCountDictTemp;
for (const auto& pair : clsCountDict) {
clsCountDictTemp.insert(std::pair<long, size_t>(pair.first, pair.second));
}
double sumdemx = 0;
for (long i = 0; i < demCount; i++) {
sumdemx= sumdemx+demX.get()[i];
}
for (long i = 0; i < demCount; i++) {
long clsid = clsArr.get()[i];
size_t Currentclscount = clsCountDictTemp[clsid];
size_t allclscount = clsCountDict[clsid];
if (clsGoalStateDict.find(clsid) == clsGoalStateDict.end()) {
continue;
}
clsGoalStateDict[clsid].get()[Currentclscount - allclscount];
clsGoalStateDict[clsid].get()[allclscount- Currentclscount].Tx = demX.get()[i];
clsGoalStateDict[clsid].get()[allclscount - Currentclscount].Ty = demY.get()[i];
clsGoalStateDict[clsid].get()[allclscount - Currentclscount].Tz = demZ.get()[i];
clsGoalStateDict[clsid].get()[allclscount - Currentclscount].TsX = slpX.get()[i];
clsGoalStateDict[clsid].get()[allclscount - Currentclscount].TsY = slpY.get()[i];
clsGoalStateDict[clsid].get()[allclscount - Currentclscount].TsZ = slpZ.get()[i];
clsGoalStateDict[clsid].get()[allclscount - Currentclscount].cls = clsArr.get()[i];
clsCountDictTemp[clsid] = clsCountDictTemp[clsid] - 1;
}
RFPCTask task;
// 参数声明
task.freqNum = this->EchoSimulationData->getPlusePoints();
task.prfNum = prfcount;
task.Rref = this->EchoSimulationData->getRefPhaseRange();
task.Rnear = this->EchoSimulationData->getNearRange();
task.Rfar = this->EchoSimulationData->getFarRange();
task.Pt = this->TaskSetting->getPt();
task.startFreq = this->EchoSimulationData->getCenterFreq() - this->EchoSimulationData->getBandwidth() / 2;
task.stepFreq = this->EchoSimulationData->getBandwidth() / (task.freqNum - 1);
task.d_echoData = (cuComplex*)mallocCUDADevice(prfcount * task.freqNum * sizeof(cuComplex), devId);
PRINT("Dev:%d ,freqnum%d , prfnum:%d ,Rref: %e ,Rnear : %e ,Rfar: %e , StartFreq: %e ,DeletFreq: %e \n",
devId,task.freqNum,task.prfNum,task.Rref,task.Rnear,task.Rfar,task.startFreq,task.stepFreq);
// 天线位置
{
std::shared_ptr<SatelliteAntPos> antplise = this->EchoSimulationData->getAntPosVelc();
std::shared_ptr<SateState> h_antlist((SateState*)mallocCUDAHost(prfcount * sizeof(SateState)), FreeCUDAHost);
for (long i = 0; i < prfcount; i++) {
h_antlist.get()[i].Px = antplise.get()[i + startprfid].Px;
h_antlist.get()[i].Py = antplise.get()[i + startprfid].Py;
h_antlist.get()[i].Pz = antplise.get()[i + startprfid].Pz;
h_antlist.get()[i].Vx = antplise.get()[i + startprfid].Vx;
h_antlist.get()[i].Vy = antplise.get()[i + startprfid].Vy;
h_antlist.get()[i].Vz = antplise.get()[i + startprfid].Vz;
h_antlist.get()[i].antDirectX = antplise.get()[i + startprfid].AntDirectX;
h_antlist.get()[i].antDirectY = antplise.get()[i + startprfid].AntDirectY;
h_antlist.get()[i].antDirectZ = antplise.get()[i + startprfid].AntDirectZ;
}
task.antlist = (SateState*)mallocCUDADevice(prfcount * sizeof(SateState), devId);
HostToDevice(h_antlist.get(), task.antlist, sizeof(SateState) * prfcount);
}
// 分块计算
for (const auto& pair : clsGoalStateDict) {
long clsid = pair.first;
size_t clscount = clsCountDict[clsid];
PRINT("Process Class ID : %d , Count: %d\n", clsid, clscount);
task.targetnum = clscount;
task.goallist = (GoalState*)mallocCUDADevice(clscount * sizeof(GoalState), devId);
HostToDevice(clsGoalStateDict[clsid].get(), task.goallist, sizeof(GoalState) * clscount);
task.sigma0_cls = clsCUDASigmaParamsDict[clsid];
ProcessRFPCTask(task,devId);
FreeCUDADevice(task.goallist);
}
this->SaveBlockSimulationEchoArr(task.d_echoData, prfcount, task.freqNum, startprfid);
FreeCUDADevice(task.d_echoData);
FreeCUDADevice(task.antlist);
//FreeCUDADevice(task.goallist);
return ErrorCode::SUCCESS;
}
ErrorCode RFPCProcessCls::SaveBlockSimulationEchoArr(cuComplex* d_echoData,size_t prfcount,size_t freqNum,long startprfid)
{
// 文件读写
omp_lock_t lock;
omp_init_lock(&lock);
omp_set_lock(&lock);
cuComplex* h_echoData = (cuComplex*)mallocCUDAHost(prfcount * freqNum * sizeof(cuComplex));
DeviceToHost(h_echoData, d_echoData, prfcount * freqNum * sizeof(cuComplex));
long prfcount_read = prfcount;
std::shared_ptr<std::complex<double>> fileEchoArr = this->EchoSimulationData->getEchoArr(startprfid, prfcount_read);
for (size_t i = 0; i < prfcount; i++) {
for (size_t j = 0; j < freqNum; j++) {
std::complex<double> temp = fileEchoArr.get()[i * freqNum + j];
fileEchoArr.get()[i * freqNum + j] = std::complex<double>(
temp.real() + h_echoData[i * freqNum + j].x,
temp.imag() + h_echoData[i * freqNum + j].y
);
}
}
qDebug() << "write echo :\t " << "prfcount:\t" << prfcount << " freqnum:\t" << freqNum;
//testOutComplexDoubleArr(QString("testoutEcho.dat"),fileEchoArr.get(), prfcount_read, freqNum);
this->EchoSimulationData->saveEchoArr(fileEchoArr, startprfid, prfcount_read);
omp_unset_lock(&lock); // 锟酵放伙拷斤拷
omp_destroy_lock(&lock); // 劫伙拷斤拷
FreeCUDAHost(h_echoData);
return ErrorCode::SUCCESS;
}