拆分回波构建，方便计算与检查回波代码，修改了回波积累结果

BaseTool/BaseConstVariable.h

@@ -2,7 +2,7 @@
 
 #ifndef BASECONSTVARIABLE_H
 #define BASECONSTVARIABLE_H
-//#define EIGEN_USE_MKL_ALL
+#define EIGEN_USE_MKL_ALL
 //#define EIGEN_NO_DEBUG 
 
  

BaseTool/EchoDataFormat.cpp

@@ -145,9 +145,10 @@ ErrorCode EchoL0Dataset::OpenOrNew(QString folder, QString filename, long PluseC
         CPLSetConfigOption("GDAL_FILENAME_IS_UTF8", "YES");
  
         GDALDriver* poDriver = GetGDALDriverManager()->GetDriverByName("ENVI");
-        std::shared_ptr<GDALDataset> poDstDS(poDriver->Create(this->GPSPointFilePath.toUtf8().constData(), 19, PluseCount, 1, GDT_Float64, NULL));
+        GDALDataset* poDstDS=(poDriver->Create(this->GPSPointFilePath.toUtf8().constData(), 19, PluseCount, 1, GDT_Float64, NULL));
-        GDALFlushCache((GDALDatasetH)poDstDS.get());
+        GDALFlushCache((GDALDatasetH)poDstDS);
-        poDstDS.reset();
+        GDALClose(poDstDS);
+        //poDstDS.reset();
         omp_unset_lock(&lock);	 //  
         omp_destroy_lock(&lock); //  
  
@@ -163,9 +164,10 @@ ErrorCode EchoL0Dataset::OpenOrNew(QString folder, QString filename, long PluseC
         CPLSetConfigOption("GDAL_FILENAME_IS_UTF8", "YES");
  
         GDALDriver* poDriver = GetGDALDriverManager()->GetDriverByName("ENVI");
-        std::shared_ptr<GDALDataset> poDstDS(poDriver->Create(this->echoDataFilePath.toUtf8().constData(), PlusePoints, PluseCount, 1, GDT_CFloat32, NULL));
+        GDALDataset* poDstDS = (poDriver->Create(this->echoDataFilePath.toUtf8().constData(), PlusePoints, PluseCount, 1, GDT_CFloat32, NULL));
-        GDALFlushCache((GDALDatasetH)poDstDS.get());
+        GDALFlushCache((GDALDatasetH)poDstDS);
-        poDstDS.reset();
+        GDALClose(poDstDS);
+        //poDstDS.reset();
         omp_unset_lock(&lock);	 //  
         omp_destroy_lock(&lock); //  
 
@@ -594,7 +596,11 @@ ErrorCode EchoL0Dataset::saveEchoArr(std::shared_ptr<std::complex<float>> echoPt
     omp_lock_t lock;
     omp_init_lock(&lock);
     omp_set_lock(&lock);
-    std::shared_ptr<GDALDataset> rasterDataset = OpenDataset(this->echoDataFilePath, GDALAccess::GA_Update);
+
+    GDALAllRegister();
+    CPLSetConfigOption("GDAL_FILENAME_IS_UTF8", "YES");
+    GDALDataset* rasterDataset = (GDALDataset*)(GDALOpen(this->echoDataFilePath.toUtf8().constData(), GDALAccess::GA_Update));
+ 
 
 
     GDALDataType	gdal_datatype = rasterDataset->GetRasterBand(1)->GetRasterDataType();
@@ -620,13 +626,14 @@ ErrorCode EchoL0Dataset::saveEchoArr(std::shared_ptr<std::complex<float>> echoPt
     else {
         if (gdal_datatype == GDT_CFloat32) {
             poBand->RasterIO(GF_Write, 0, startPRF, width, PRFLen,  echoPtr.get(), width, PRFLen, GDT_CFloat32, 0, 0);
-            GDALFlushCache((GDALDatasetH)rasterDataset.get());
+            GDALFlushCache((GDALDatasetH)rasterDataset);
         }
         else {
             qDebug() << QString::fromStdString(errorCode2errInfo(ErrorCode::ECHO_L0DATA_ECHOFILEFORMATERROR));
         }
     }
-    rasterDataset.reset();
+    GDALClose(rasterDataset);
+    rasterDataset = nullptr;
     omp_unset_lock(&lock);	 //  
     omp_destroy_lock(&lock); //  
     return ErrorCode::SUCCESS;

BaseTool/EchoDataFormat.h

@@ -182,10 +182,6 @@ public:
 	ErrorCode saveAntPos(std::shared_ptr<double> ptr); // 注意这个方法很危险，请写入前检查数据是否正确
 	ErrorCode saveEchoArr(std::shared_ptr<std::complex<float>> echoPtr, long startPRF, long PRFLen);
 
-
-
-
-
 };
 
 

GPUTool/GPURTPC.cu

@@ -112,16 +112,13 @@ __device__ float GPU_BillerInterpAntPattern(float* antpattern,
 
 __device__ cuComplex  GPU_calculationEcho(float sigma0, float TransAnt, float ReciveAnt,
 	float localangle, float R, float slopeangle, float Pt, float lamda) {
-	float r = R;
 	float amp = Pt * TransAnt * ReciveAnt;
 	amp = amp * sigma0;
-	amp = amp / (powf(4 * LAMP_CUDA_PI, 2) * powf(r, 4)); // ·´ÉäÇ¿¶È
+	amp = amp / (powf(4 * LAMP_CUDA_PI, 2) * powf(R, 4)); // 反射强度
-	float phi = (-4 * LAMP_CUDA_PI / lamda) * r;
+	float phi = (-4 * LAMP_CUDA_PI / lamda) * R;
 	cuComplex echophi = make_cuComplex(0, phi);
 	cuComplex echophiexp = cuCexpf(echophi);
-	cuComplex echo;
+	cuComplex echo=make_cuComplex(echophiexp.x * amp, echophiexp.y * amp);
-	echo.x = echophiexp.x * amp;
-	echo.y = echophiexp.y * amp;
 	return echo;
 }
 
@@ -295,14 +292,12 @@ __global__ void CUDA_calculationEcho(float* sigma0, float* TransAnt, float* Reci
 
 		float timeR = 2 * (r - nearRange) / LIGHTSPEED * Fs;
 		long timeID = floorf(timeR);
-		if (timeID < 0 || timeID >= FreqIDmax) {
+		//if (timeID < 0 || timeID >= FreqIDmax) {
-			timeID = 0;
+		//	timeID = 0;
-			amp = 0;
+		//	amp = 0;
-		}
+		//}
 
-		cuComplex echo;
+		cuComplex echo = make_cuComplex(echophiexp.x , echophiexp.y);
-		echo.x = echophiexp.x * amp;
-		echo.y = echophiexp.y * amp;
 		echoArr[idx] = echo;
 		FreqID[idx] = timeID;
 	}
@@ -317,8 +312,6 @@ __global__ void CUDA_AntPatternInterpGain(float* anttheta, float* antphi, float*
 
 		float temptheta = anttheta[idx];
 		float tempphi = antphi[idx];
-		
-
 		float antPatternGain =  GPU_BillerInterpAntPattern(antpattern,
 			starttheta, startphi, dtheta, dphi, thetapoints, phipoints,
 			temptheta, tempphi) ;
@@ -361,6 +354,56 @@ __global__ void CUDA_InterpSigma(
 }
 
 
+__global__ void CUDA_CalculationEchoAmp(float* sigma0, float* TransAnt, float* ReciveAnt, float* R, 
+	float Pt, 
+	float* ampArr, long len) {
+	long  idx = blockIdx.x * blockDim.x + threadIdx.x;
+	if (idx < len) {
+		float r = R[idx];
+		float amptemp = Pt * TransAnt[idx] * ReciveAnt[idx] * sigma0[idx];
+		amptemp = amptemp / (powf(4 * LAMP_CUDA_PI, 2) * powf(r, 4)); // 反射强度
+		ampArr[idx] = amptemp;
+	}
+}
+
+
+__global__ void CUDA_CalculationEchoPhase(float* R, float lamda, float* phaseArr, long len) {
+	long  idx = blockIdx.x * blockDim.x + threadIdx.x;
+	if (idx < len) {
+		float r = R[idx];
+		// 处理相位
+		float phi = (-4 * LAMP_CUDA_PI / lamda) * r;
+		phaseArr[idx] = phi;
+	}
+}
+
+
+__global__ void CUDA_CombinationEchoAmpAndPhase(float* R,
+	float* echoAmp,float* echoPhase,
+	float nearRange, float Fs, long plusepoints,
+	cuComplex* echo, long* FreqID, long len) {
+	long  idx = blockIdx.x * blockDim.x + threadIdx.x;
+	if (idx < len) {
+		float r = R[idx];
+		float phase = echoPhase[idx];
+		float amp = echoAmp[idx];
+		cuComplex echophi = make_cuComplex(0, phase);
+		cuComplex echophiexp = cuCexpf(echophi);
+
+		float timeR = 2 * (r - nearRange) / LIGHTSPEED * Fs;
+		long timeID = floorf(timeR);
+		if (timeID < 0 || timeID >= plusepoints) {
+			timeID = 0;
+			amp = 0;
+		}
+		cuComplex echotemp = make_cuComplex(echophiexp.x*amp, echophiexp.y*amp);
+		echo[idx] = echotemp;
+		FreqID[idx] = timeID;
+
+	}
+}
+
+
 
 extern "C" void SatelliteAntDirectNormal(float* RstX, float* RstY, float* RstZ,
 	float  antXaxisX, float  antXaxisY, float  antXaxisZ,
@@ -437,6 +480,68 @@ extern "C" void calculationEcho(float* sigma0, float* TransAnt, float* ReciveAnt
 	cudaDeviceSynchronize();
 }
 
+void CUDACalculationEchoAmp(float* sigma0, float* TransAnt, float* ReciveAnt, float* R, float Pt, float* ampArr, long len)
+{
+
+	int blockSize = 256; // 每个块的线程数
+	int numBlocks = (len + blockSize - 1) / blockSize; // 根据 pixelcount 计算网格大小
+	// 调用 CUDA 核函数
+	CUDA_CalculationEchoAmp << <numBlocks, blockSize >> > (
+		sigma0, TransAnt, ReciveAnt, R,	Pt, ampArr, len);
+#ifdef __CUDADEBUG__
+	cudaError_t err = cudaGetLastError();
+	if (err != cudaSuccess) {
+		printf("CUDA_RTPC_SiglePRF CUDA Error: %s\n", cudaGetErrorString(err));
+		// Possibly: exit(-1) if program cannot continue....
+	}
+#endif // __CUDADEBUG__
+	cudaDeviceSynchronize();
+
+}
+
+void CUDACalculationEchoPhase(float* R, float lamda, float* phaseArr, long len)
+{
+	int blockSize = 256; // 每个块的线程数
+	int numBlocks = (len + blockSize - 1) / blockSize; // 根据 pixelcount 计算网格大小
+	// 调用 CUDA 核函数
+	CUDA_CalculationEchoPhase << <numBlocks, blockSize >> > (
+		R, lamda, phaseArr, len);
+#ifdef __CUDADEBUG__
+	cudaError_t err = cudaGetLastError();
+	if (err != cudaSuccess) {
+		printf("CUDA_RTPC_SiglePRF CUDA Error: %s\n", cudaGetErrorString(err));
+		// Possibly: exit(-1) if program cannot continue....
+	}
+#endif // __CUDADEBUG__
+	cudaDeviceSynchronize();
+
+
+}
+
+void CUDACombinationEchoAmpAndPhase(float* R,
+	float* echoAmp, float* echoPhase,
+	float nearRange, float Fs, long plusepoints, cuComplex* echo, long* FreqID, long len)
+{
+
+	int blockSize = 256; // 每个块的线程数
+	int numBlocks = (len + blockSize - 1) / blockSize; // 根据 pixelcount 计算网格大小
+	// 调用 CUDA 核函数
+	CUDA_CombinationEchoAmpAndPhase << <numBlocks, blockSize >> > (
+		R, 
+		echoAmp, echoPhase,
+		nearRange, Fs, plusepoints, echo, FreqID, len
+		);
+#ifdef __CUDADEBUG__
+	cudaError_t err = cudaGetLastError();
+	if (err != cudaSuccess) {
+		printf("CUDA_RTPC_SiglePRF CUDA Error: %s\n", cudaGetErrorString(err));
+		// Possibly: exit(-1) if program cannot continue....
+	}
+#endif // __CUDADEBUG__
+	cudaDeviceSynchronize();
+
+}
+
 
 extern "C" void CUDAInterpSigma(
 	long* demcls,float* sigmaAmp, float* localanglearr,long len,

GPUTool/GPURTPC.cuh

@@ -36,6 +36,28 @@ extern "C" void calculationEcho(float* sigma0, float* TransAnt, float* ReciveAnt
 	cuComplex* echoAmp, long* FreqID,
 	long len);
 
+extern "C" void CUDACalculationEchoAmp(
+	float* sigma0, 
+	float* TransAnt, float* ReciveAnt,
+	float* R,
+	float Pt,
+	float* ampArr,
+	long len
+);
+
+extern "C" void CUDACalculationEchoPhase(
+	float* R, float lamda,
+	float* phaseArr,
+	long len
+);
+
+
+extern "C" void CUDACombinationEchoAmpAndPhase(float* R, float* echoAmp, float* echoPhase,
+	float nearRange, float Fs,long plusepoints,
+	cuComplex* echo,long* FreqID,
+	long len
+);
+
 
 extern "C" void CUDAInterpSigma(
 	long* demcls, float* sigmaAmp, float* localanglearr, long len,

GPUTool/GPUTBPImage.cu

@@ -19,7 +19,7 @@
 
 __global__ void CUDA_TBPImage(
 	float* antPx, float* antPy, float* antPz,
-	float* imgx, float* imgy, float* imgz,
+	float* imgx, float* imgy, float* imgz, float* RArr,
 	cuComplex* echoArr, cuComplex* imgArr,
 	float freq, float fs, float Rnear, float Rfar,
 	long rowcount, long colcount,
@@ -31,11 +31,16 @@ __global__ void CUDA_TBPImage(
 		float R = sqrtf(powf(antPx[prfid] - imgx[idx], 2) + powf(antPy[prfid] - imgy[idx], 2) + powf(antPz[prfid] - imgz[idx], 2));
 		float Ridf = 2 * (R - Rnear) / LIGHTSPEED * fs;
 		long Rid = floorf(Ridf);
+		RArr[idx] = R;
 		if(Rid <0|| Rid >= freqcount){}
 		else {
 			float factorj = freq * 4 * PI / LIGHTSPEED;
-			cuComplex Rphi =cuCexpf(make_cuComplex(0, factorj * R));// ĐŁŐýĎî
+			cuComplex Rfactorj = make_cuComplex(0, factorj * R);
+			cuComplex Rphi =cuCexpf(Rfactorj);// ĐŁŐýĎî
 			imgArr[idx] = cuCaddf(imgArr[idx], cuCmulf(echoArr[Rid] , Rphi));// ½ÃÕý
+			//printf("R=%f;Rid=%d;factorj=%f;Rfactorj=complex(%f,%f);Rphi=complex(%f,%f);\n", R, Rid, factorj, 
+			//	Rfactorj.x, Rfactorj.y,
+			//	Rphi.x, Rphi.y);
 		}
 	}
 }
@@ -43,6 +48,7 @@ __global__ void CUDA_TBPImage(
 
 extern "C" void CUDATBPImage(float* antPx, float* antPy, float* antPz,
 	float* imgx, float* imgy, float* imgz, 
+	float* R,
 	cuComplex* echoArr, cuComplex* imgArr,
 	float freq, float fs, float Rnear, float Rfar,
 	long rowcount, long colcount, 
@@ -54,8 +60,9 @@ extern "C" void CUDATBPImage(float* antPx, float* antPy, float* antPz,
 	// µ÷Óà CUDA ºËº¯Êý CUDA_RTPC_Kernel
 
 	CUDA_TBPImage << <numBlocks, blockSize >> > (
-		 antPx,   antPy,   antPz,
+			 antPx,   antPy,   antPz,
 		  imgx,  imgy,   imgz,
+			R,
 		  echoArr,   imgArr,
 		  freq,   fs,   Rnear,   Rfar,
 		  rowcount,   colcount,

GPUTool/GPUTBPImage.cuh

@@ -13,7 +13,7 @@
 
 extern __global__ void CUDA_TBPImage(
 	float* antPx, float* antPy, float* antPz,
-	float* imgx, float* imgy, float* imgz,
+	float* imgx, float* imgy, float* imgz, float* R,
 	cuComplex* echoArr, cuComplex* imgArr,
 	float freq, float fs, float Rnear, float Rfar,
 	long rowcount, long colcount,
@@ -27,6 +27,7 @@ extern "C" void CUDATBPImage(
 	float* imgx,
 	float* imgy,
 	float* imgz,
+	float* R,
 	cuComplex* echoArr,
 	cuComplex* imgArr,
 	float freq, float fs, float Rnear, float Rfar,

GPUTool/GPUTool.cu

@@ -90,6 +90,30 @@ __global__ void CUDA_cosAngle_VA_AB(float* Ax, float* Ay, float* Az, float* Bx,
 	}
 }
 
+__global__ void CUDA_GridPoint_Linear_Interp1(float* v, float* q, float* qv, long xlen, long qlen)
+{
+	long  idx = blockIdx.x * blockDim.x + threadIdx.x;
+	if (idx < qlen) {
+		float qx = q[idx];
+		// 检索循环
+		if (qx < 0 || qx > xlen - 1) {}
+		else {
+			long x1 = floor(qx);
+			long x2 = ceil(qx);
+
+			if (x1 >= 0 && x2 < xlen) {
+				float y1 = v[x1];
+				float y2 = v[x2];
+				float y = y1 + (y2 - y1) * (qx - x1) / (x2 - x1);
+				qv[idx] = y;
+			}
+			else {
+
+			}
+		}
+	}
+}
+
 //´íÎóÌáʾ
 extern "C"  void checkCudaError(cudaError_t err, const char* msg) {
 	if (err != cudaSuccess) {
@@ -194,7 +218,7 @@ extern "C" void CUDAdistanceAB(float* Ax, float* Ay, float* Az, float* Bx, float
 #ifdef __CUDADEBUG__
 	cudaError_t err = cudaGetLastError();
 	if (err != cudaSuccess) {
-		printf("CUDA_RTPC_SiglePRF CUDA Error: %s\n", cudaGetErrorString(err));
+		printf("CUDAdistanceAB CUDA Error: %s\n", cudaGetErrorString(err));
 		// Possibly: exit(-1) if program cannot continue....
 	}
 #endif // __CUDADEBUG__
@@ -211,7 +235,7 @@ extern "C" void CUDABdistanceAs(float* Ax, float* Ay, float* Az, float Bx, float
 #ifdef __CUDADEBUG__
 	cudaError_t err = cudaGetLastError();
 	if (err != cudaSuccess) {
-		printf("CUDA_RTPC_SiglePRF CUDA Error: %s\n", cudaGetErrorString(err));
+		printf("CUDABdistanceAs CUDA Error: %s\n", cudaGetErrorString(err));
 		// Possibly: exit(-1) if program cannot continue....
 	}
 #endif // __CUDADEBUG__
@@ -228,7 +252,7 @@ extern "C" void CUDAmake_VectorA_B(float sX, float sY, float sZ, float* tX, floa
 #ifdef __CUDADEBUG__
 	cudaError_t err = cudaGetLastError();
 	if (err != cudaSuccess) {
-		printf("CUDA_RTPC_SiglePRF CUDA Error: %s\n", cudaGetErrorString(err));
+		printf("CUDAmake_VectorA_B CUDA Error: %s\n", cudaGetErrorString(err));
 		// Possibly: exit(-1) if program cannot continue....
 	}
 #endif // __CUDADEBUG__
@@ -245,7 +269,7 @@ extern "C" void CUDANorm_Vector(float* Vx, float* Vy, float* Vz, float* R, long
 #ifdef __CUDADEBUG__
 	cudaError_t err = cudaGetLastError();
 	if (err != cudaSuccess) {
-		printf("CUDA_RTPC_SiglePRF CUDA Error: %s\n", cudaGetErrorString(err));
+		printf("CUDANorm_Vector CUDA Error: %s\n", cudaGetErrorString(err));
 		// Possibly: exit(-1) if program cannot continue....
 	}
 #endif // __CUDADEBUG__
@@ -261,13 +285,32 @@ extern "C" void CUDAcosAngle_VA_AB(float* Ax, float* Ay, float* Az, float* Bx, f
 #ifdef __CUDADEBUG__
 	cudaError_t err = cudaGetLastError();
 	if (err != cudaSuccess) {
-		printf("CUDA_RTPC_SiglePRF CUDA Error: %s\n", cudaGetErrorString(err));
+		printf("CUDAcosAngle_VA_AB CUDA Error: %s\n", cudaGetErrorString(err));
 		// Possibly: exit(-1) if program cannot continue....
 	}
 #endif // __CUDADEBUG__
 	cudaDeviceSynchronize();
 }
 
+extern "C" void  CUDAGridPointLinearInterp1(float* v, float* q, float* qv, long xlen, long qlen)
+{
+
+	int blockSize = 256; // 每个块的线程数
+	int numBlocks = (qlen + blockSize - 1) / blockSize; // 根据 pixelcount 计算网格大小
+	// 调用 CUDA 核函数
+	CUDA_GridPoint_Linear_Interp1 << <numBlocks, blockSize >> > ( v, q,qv, xlen, qlen);
+
+#ifdef __CUDADEBUG__
+	cudaError_t err = cudaGetLastError();
+	if (err != cudaSuccess) {
+		printf("CUDALinearInterp1 CUDA Error: %s\n", cudaGetErrorString(err));
+		// Possibly: exit(-1) if program cannot continue....
+	}
+#endif // __CUDADEBUG__
+	cudaDeviceSynchronize();
+
+}
+
 #endif
 
 

GPUTool/GPUTool.cuh

@@ -49,7 +49,7 @@ extern __global__ void CUDA_B_DistanceA(float* Ax, float* Ay, float* Az, float B
 extern __global__ void CUDA_make_VectorA_B(float sX, float sY, float sZ, float* tX, float* tY, float* tZ, float* RstX, float* RstY, float* RstZ, long len);
 extern __global__ void CUDA_Norm_Vector(float* Vx, float* Vy, float* Vz, float* R, long len);
 extern __global__ void CUDA_cosAngle_VA_AB(float* Ax, float* Ay, float* Az, float* Bx, float* By, float* Bz, float* anglecos, long len);
-
+extern __global__ void CUDA_GridPoint_Linear_Interp1(float* v, float* q, float* qv, long xlen, long qlen);
 
 // 误差处理函数
 extern "C"  void checkCudaError(cudaError_t err, const char* msg);
@@ -63,12 +63,19 @@ extern "C" void HostToDevice(void* hostptr, void* deviceptr, long memsize);//GPU
 extern "C" void DeviceToHost(void* hostptr, void* deviceptr, long memsize);//GPU 内存数据转移  GPU -> 设备
 
 
-// 基础运算函数
+// 矢量基础运算函数
 extern "C" void CUDAdistanceAB(float* Ax, float* Ay, float* Az, float* Bx, float* By, float* Bz, float* R, long member);
 extern "C" void CUDABdistanceAs(float* Ax, float* Ay, float* Az, float Bx, float By, float Bz, float* R, long member);
 extern "C" void CUDAmake_VectorA_B(float sX, float sY, float sZ, float* tX, float* tY, float* tZ, float* RstX, float* RstY, float* RstZ, long member);
 extern "C" void CUDANorm_Vector(float* Vx, float* Vy, float* Vz, float* R, long member);
 extern "C" void CUDAcosAngle_VA_AB(float* Ax, float* Ay, float* Az, float* Bx, float* By, float* Bz, float* anglecos, long len);
 
+// 常见插值算法
+
+extern "C" void CUDAGridPointLinearInterp1(float* v, float* q, float* qv,long xlen, long qlen);
+
+
+
+
 #endif
 #endif

Imageshow/ImageShowDialogClass.cpp

@@ -207,6 +207,35 @@ void ImageShowDialogClass::load_double_MatrixX_data(Eigen::MatrixXd X, Eigen::Ma
 	ui->m_plot->replot();
 }
 
+void ImageShowDialogClass::load_double_data(float* data, long rows, long cols, QString name)
+{
+	Eigen::MatrixXd dataArr = Eigen::MatrixXd::Zero(rows, cols);
+	for (long i = 0; i < rows; i++) {
+		for (long j = 0; j < cols; j++) {
+			dataArr(i, j) = data[i * cols + j];
+		}
+	}
+	this->load_double_MatrixX_data(dataArr, name);
+}
+
+void ImageShowDialogClass::load_double_data(float* Xs, float* Ys, float* data, long rows, long cols, QString name)
+{
+
+	Eigen::MatrixXd dataArr = Eigen::MatrixXd::Zero(rows, cols);
+	Eigen::MatrixXd dataX = Eigen::MatrixXd::Zero(rows, cols);
+	Eigen::MatrixXd dataY = Eigen::MatrixXd::Zero(rows, cols);
+	for (long i = 0; i < rows; i++) {
+		for (long j = 0; j < cols; j++) {
+			dataArr(i, j) = data[i * cols + j];
+			dataX(i, j) = Xs[i * cols + j];
+			dataY(i, j) = Ys[i * cols + j];
+		}
+	}
+	this->load_double_MatrixX_data(dataX,dataY,dataArr, name);
+
+
+}
+
 void ImageShowDialogClass::remove_Data(QString name)
 {
 	for(size_t i=0;i<this->m_plot->graphCount();i++){

Imageshow/ImageShowDialogClass.h

@@ -132,6 +132,9 @@ public:
 	void Scatter(QVector<double> xs, QVector<double> ys);
 	void load_double_MatrixX_data(Eigen::MatrixXd data, QString name);
 	void load_double_MatrixX_data(Eigen::MatrixXd X,Eigen::MatrixXd Y,Eigen::MatrixXd data, QString name);
+	void load_double_data(float* data, long rows, long cols, QString name);
+	void load_double_data(float* Xs,float* Ys,float* data, long rows, long cols, QString name);
+
 	void remove_Data( QString name);
 
 	LAMPDATASHOWCLASS getGraphClass(size_t i = 0);

SimulationSAR/QImageSARRTPC.ui

@@ -77,7 +77,7 @@
           </size>
          </property>
          <property name="text">
-          <string>D:/Programme/vs2022/RasterMergeTest/TestData/ant/ant_model_setting_Horn_conical1_FarField-receive.csv</string>
+          <string>D:/Programme/vs2022/RasterMergeTest/simulationData/ant/ant_model_setting_Horn_conical1_FarField-receive.csv</string>
          </property>
         </widget>
        </item>
@@ -90,7 +90,7 @@
           </size>
          </property>
          <property name="text">
-          <string>D:/Programme/vs2022/RasterMergeTest/TestData/landcover_aligned2.dat</string>
+          <string>D:/Programme/vs2022/RasterMergeTest/simulationData/landcover_aligned2.dat</string>
          </property>
         </widget>
        </item>
@@ -116,7 +116,7 @@
           </size>
          </property>
          <property name="text">
-          <string>D:/Programme/vs2022/RasterMergeTest/TestData/GF3_Simulation_Setting.xml</string>
+          <string>D:/Programme/vs2022/RasterMergeTest/simulationData/GF3_Simulation_Setting.xml</string>
          </property>
         </widget>
        </item>
@@ -142,7 +142,7 @@
           </size>
          </property>
          <property name="text">
-          <string>D:/Programme/vs2022/RasterMergeTest/TestData/ant/ant_model_setting_Horn_conical1_FarField-trans.csv</string>
+          <string>D:/Programme/vs2022/RasterMergeTest/simulationData/ant/ant_model_setting_Horn_conical1_FarField-trans.csv</string>
          </property>
         </widget>
        </item>
@@ -194,7 +194,7 @@
           </size>
          </property>
          <property name="text">
-          <string>D:/Programme/vs2022/RasterMergeTest/TestData/115E39N_COP30_clip.tif</string>
+          <string>D:/Programme/vs2022/RasterMergeTest/simulationData/115E39N_COP30_clip.tif</string>
          </property>
         </widget>
        </item>
@@ -272,7 +272,7 @@
           </size>
          </property>
          <property name="text">
-          <string>D:/Programme/vs2022/RasterMergeTest/TestData/GF3_Simulation_GPSNode.xml</string>
+          <string>D:/Programme/vs2022/RasterMergeTest/simulationData/GF3_Simulation_GPSNode.xml</string>
          </property>
         </widget>
        </item>

SimulationSAR/RTPCProcessCls.cpp

@@ -132,6 +132,13 @@ ErrorCode RTPCProcessCls::Process(long num_thread)
 	else {}
 	qDebug() << "RTPCMainProcess";
 
+	stateCode = this->InitEchoMaskArray();
+	if (stateCode != ErrorCode::SUCCESS) {
+		return stateCode;
+	}
+	else {}
+	qDebug() << "InitEchoMaskArray";
+
 	//stateCode = this->RTPCMainProcess(num_thread);
 	stateCode = this->RTPCMainProcess_GPU( );
 
@@ -188,9 +195,9 @@ ErrorCode RTPCProcessCls::InitParams()
 ErrorCode RTPCProcessCls::DEMPreprocess()
 {
 
-	this->demxyzPath = QDir(tmpfolderPath).filePath("demxyz.tif");
+	this->demxyzPath = QDir(tmpfolderPath).filePath("demxyz.bin");
 	gdalImage demds(this->DEMTiffPath);
-	gdalImage demxyz = CreategdalImage(demxyzPath, demds.height, demds.width, 3, demds.gt, demds.projection, true, true);// X,Y,Z
+	gdalImage demxyz = CreategdalImage(demxyzPath, demds.height, demds.width, 3, demds.gt, demds.projection, true, true,true);// X,Y,Z
 
 	// 分块计算并转换为XYZ
 
@@ -234,11 +241,11 @@ ErrorCode RTPCProcessCls::DEMPreprocess()
 
 	// 计算坡向角
 
-	this->demsloperPath = QDir(tmpfolderPath).filePath("demsloper.tif");
+	this->demsloperPath = QDir(tmpfolderPath).filePath("demsloper.bin");
-	this->demmaskPath = QDir(tmpfolderPath).filePath("demmask.tif");
+	this->demmaskPath = QDir(tmpfolderPath).filePath("demmask.bin");
 
-	gdalImage demsloperxyz = CreategdalImage(this->demsloperPath, demds.height, demds.width, 4, demds.gt, demds.projection, true, true);// X,Y,Z,cosangle
+	gdalImage demsloperxyz = CreategdalImage(this->demsloperPath, demds.height, demds.width, 4, demds.gt, demds.projection, true, 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
+	gdalImage demmask = CreategdalImage(this->demmaskPath, demxyz.height, demxyz.width, 1, demxyz.gt, demxyz.projection, true, true,true);// X,Y,Z
 
 	line_invert = 1000;
 	long start_ids = 0;
@@ -300,6 +307,34 @@ ErrorCode RTPCProcessCls::DEMPreprocess()
 	return ErrorCode::SUCCESS;
 }
 
+ErrorCode RTPCProcessCls::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;
+}
+
 
 ErrorCode RTPCProcessCls::RTPCMainProcess_GPU( )
 {
@@ -318,61 +353,11 @@ ErrorCode RTPCProcessCls::RTPCMainProcess_GPU( )
 	imageStarttime = this->TaskSetting->getSARImageStartTime();
 	//std::vector<SatelliteOribtNode> sateOirbtNodes(this->PluseCount);
 
-	std::shared_ptr<SatelliteOribtNode[]> sateOirbtNodes(new SatelliteOribtNode[this->PluseCount], delArrPtr);
+	std::shared_ptr<SatelliteOribtNode[]> sateOirbtNodes = this->getSatelliteOribtNodes(prf_time, dt, antflag, imageStarttime); // 获取天线坐标
-	{ // 姿态计算不同
-		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 !!!";
-	}
 
 	// 回波
 	long echoIdx = 0;
-	double NearRange = this->EchoSimulationData->getNearRange(); // 近斜据
+	double NearRange = this->EchoSimulationData->getNearRange(); // 近斜距
 	double FarRange = this->EchoSimulationData->getFarRange();
 
 	double TimgNearRange = 2 * NearRange / LIGHTSPEED;
@@ -393,6 +378,9 @@ ErrorCode RTPCProcessCls::RTPCMainProcess_GPU( )
 	this->EchoSimulationData->initEchoArr(std::complex<float>(0, 0));
 	POLARTYPEENUM polartype = this->TaskSetting->getPolarType();
  
+	gdalImage echoMaskImg(this->OutEchoMaskPath);
+
+
 #ifndef __CUDANVCC___ 
 	QMessageBox::information(this, u8"程序提示", u8"请确定安装了CUDA库");
 #else
@@ -544,13 +532,9 @@ ErrorCode RTPCProcessCls::RTPCMainProcess_GPU( )
 		std::cout << "";
 	}
 
-
-
 	HostToDevice(h_clsSigmaParam, d_clsSigmaParam, sizeof(CUDASigmaParam) * clamapid);
 	
-
+	// 地面 XYZ
-
-
 	Eigen::MatrixXd dem_x = demxyz.getData(0, 0, blokline, demxyz.width, 1);  // 地面坐标
 	long tempDemRows = dem_x.rows();
 	long tempDemCols = dem_x.cols();
@@ -654,10 +638,10 @@ ErrorCode RTPCProcessCls::RTPCMainProcess_GPU( )
 	d_ReciveAnt= (float*)mallocCUDADevice(  sizeof(float) * blokline * tempDemCols);
 
 	// 回波
-	cuComplex* h_echoAmp;
+	cuComplex* h_echo;
-	cuComplex* d_echoAmp;
+	cuComplex* d_echo;
-	h_echoAmp=(cuComplex*)mallocCUDAHost(sizeof(cuComplex) * blokline * tempDemCols);
+	h_echo=(cuComplex*)mallocCUDAHost(sizeof(cuComplex) * blokline * tempDemCols);
-	d_echoAmp=(cuComplex*)mallocCUDADevice( sizeof(cuComplex) * blokline * tempDemCols); //19
+	d_echo=(cuComplex*)mallocCUDADevice( sizeof(cuComplex) * blokline * tempDemCols); //19
 
 	long* h_FreqID;
 	long* d_FreqID;
@@ -672,6 +656,17 @@ ErrorCode RTPCProcessCls::RTPCMainProcess_GPU( )
 	float* h_amp=(float*)mallocCUDAHost(sizeof(float) * blokline * tempDemCols);
 	float* d_amp=(float*)mallocCUDADevice(sizeof(float) * blokline * tempDemCols);
 
+	float* h_echoAmp = (float*)mallocCUDAHost(sizeof(float) * blokline * tempDemCols);
+	float* d_echoAmp = (float*)mallocCUDADevice(sizeof(float) * blokline * tempDemCols);
+
+
+	float* h_echoPhase = (float*)mallocCUDAHost(sizeof(float) * blokline * tempDemCols);
+	float* d_echoPhase = (float*)mallocCUDADevice(sizeof(float) * blokline * tempDemCols);
+
+	//float* h_R = (float*)mallocCUDAHost(sizeof(float) * blokline * tempDemCols);
+	//float* d_R = (float*)mallocCUDADevice(sizeof(float) * blokline * tempDemCols);
+
+
 	for (startline = 0; startline < demRow; startline = startline + blokline) {
 		long newblokline = blokline;
 		if ((startline + blokline) >= demRow) {
@@ -706,11 +701,14 @@ ErrorCode RTPCProcessCls::RTPCMainProcess_GPU( )
 			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_echo);						FreeCUDADevice(d_echo);//19
 			FreeCUDAHost(h_FreqID);						FreeCUDADevice(d_FreqID);//20
 			FreeCUDAHost(h_demcls);						FreeCUDADevice(d_demcls);
 			FreeCUDAHost(h_amp);						FreeCUDADevice(d_amp);
+			FreeCUDAHost(h_echoAmp);					FreeCUDADevice(d_echoAmp);
+			FreeCUDAHost(h_echoPhase);					FreeCUDADevice(d_echoPhase);
 
+		 
 			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);
@@ -728,10 +726,15 @@ ErrorCode RTPCProcessCls::RTPCMainProcess_GPU( )
 			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_echo = (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);
+			h_R = (float*)mallocCUDAHost(sizeof(float) * newblokline * tempDemCols);
+			h_echoAmp = (float*)mallocCUDAHost(sizeof(float) * newblokline * tempDemCols);
+			h_echoPhase = (float*)mallocCUDAHost(sizeof(float) * newblokline * tempDemCols);
+			
+
 
 			d_dem_x=(float*)mallocCUDADevice(sizeof(float) * newblokline * tempDemCols);
 			d_dem_y=(float*)mallocCUDADevice(sizeof(float) * newblokline * tempDemCols);
@@ -750,10 +753,13 @@ ErrorCode RTPCProcessCls::RTPCMainProcess_GPU( )
 			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_echo=(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);
+			d_R = (float*)mallocCUDADevice(sizeof(float) * newblokline * tempDemCols);
+			d_echoAmp = (float*)mallocCUDADevice(sizeof(float) * newblokline * tempDemCols);
+			d_echoPhase = (float*)mallocCUDADevice(sizeof(float) * newblokline * tempDemCols);
 		}
 //#   pragma omp parallel for
 		for (long i = 0; i < newblokline; i++) {
@@ -767,6 +773,8 @@ ErrorCode RTPCProcessCls::RTPCMainProcess_GPU( )
 				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;
+				h_echoAmp[i * demxyz.width + j] = 0.0f;
+				h_echoPhase[i * demxyz.width + j] = 0.0f;
 			}
 		}
 
@@ -781,6 +789,8 @@ ErrorCode RTPCProcessCls::RTPCMainProcess_GPU( )
 		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);
+		HostToDevice((void*)h_echoAmp, (void*)d_echoAmp, sizeof(float) * newblokline* tempDemCols);
+		HostToDevice((void*)h_echoPhase, (void*)d_echoPhase, sizeof(float) * newblokline* tempDemCols);
 
 
 
@@ -801,7 +811,7 @@ ErrorCode RTPCProcessCls::RTPCMainProcess_GPU( )
 		for (startprfid = 0; startprfid < pluseCount; startprfid = startprfid + echoblockline) {
 			long templine = startprfid + echoblockline < PluseCount ? echoblockline : PluseCount - startprfid;
 			std::shared_ptr<std::complex<float>> echotemp = this->EchoSimulationData->getEchoArr(startprfid, templine);
-
+			Eigen::MatrixXd echoMasktemp = echoMaskImg.getData(startprfid,0, templine,echoMaskImg.width,1);
 			for (long tempprfid = 0; tempprfid < templine; tempprfid++) {
 				{// 计算
 					long prfid = tempprfid + startprfid;
@@ -839,13 +849,20 @@ ErrorCode RTPCProcessCls::RTPCMainProcess_GPU( )
 						d_dem_theta, d_dem_phi, pixelcount);// 计算角度
 
 #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);
+					CUDA_Norm_Vector(d_RstX, d_RstY, d_RstZ, d_R, pixelcount);
-					//testOutAmpArr("h_RstY.bin", h_RstY, newblokline, tempDemCols);
+
-					//testOutAmpArr("h_RstZ.bin", h_RstZ, newblokline, tempDemCols);
+
+					DeviceToHost(h_RstX, d_RstX, sizeof(float)* pixelcount);
+					DeviceToHost(h_RstY, d_RstY, sizeof(float)* pixelcount);
+					DeviceToHost(h_RstZ, d_RstZ, sizeof(float)* pixelcount);
+					DeviceToHost(h_R, d_R, sizeof(float)* pixelcount);
+
+					testOutAmpArr(QString("h_RstX_%1.bin").arg(prfid), h_RstX, newblokline, tempDemCols);
+					testOutAmpArr(QString("h_RstY_%1.bin").arg(prfid), h_RstY, newblokline, tempDemCols);
+					testOutAmpArr(QString("h_RstZ_%1.bin").arg(prfid), h_RstZ, newblokline, tempDemCols);
+					testOutAmpArr(QString("h_R_%1.bin").arg(prfid), h_R, 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);
@@ -864,17 +881,25 @@ ErrorCode RTPCProcessCls::RTPCMainProcess_GPU( )
 					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);
+					// 计算强度
+					CUDACalculationEchoAmp(d_amp, d_TransAnt, d_ReciveAnt, d_R, Pt, d_echoAmp, pixelcount);
+					// 计算相位
+					CUDACalculationEchoPhase(d_R,lamda, d_echoPhase, pixelcount);
+					// 计算回波
+					CUDACombinationEchoAmpAndPhase(d_R, d_echoAmp, d_echoPhase,NearRange,Fs,PlusePoint,d_echo,d_FreqID,pixelcount);
 #ifdef __PRFDEBUG__	
-					DeviceToHost(h_amp, d_amp, sizeof(float)* pixelcount);
+					DeviceToHost(h_amp, d_amp, sizeof(float) * pixelcount);
-					testOutAmpArr(QString("amp_%1.bin").arg(prfid), h_amp,newblokline, tempDemCols);
+					DeviceToHost(h_echoAmp, d_echoAmp, sizeof(float) * pixelcount);
-					
+					DeviceToHost(h_echoPhase, d_echoPhase, sizeof(float) * pixelcount);
 
+					testOutAmpArr(QString("amp_%1.bin").arg(prfid), h_amp, newblokline, tempDemCols);
+					testOutAmpArr(QString("echoAmp_%1.bin").arg(prfid), h_echoAmp, newblokline, tempDemCols);
+					testOutAmpArr(QString("echoPhase_%1.bin").arg(prfid), h_echoPhase, newblokline, tempDemCols);
 #endif // __PRFDEBUG__ 
 
-					// 计算回波
+					DeviceToHost(h_echo, d_echo, sizeof(cuComplex) * pixelcount);
-					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);
 
@@ -899,8 +924,12 @@ ErrorCode RTPCProcessCls::RTPCMainProcess_GPU( )
 #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);
+						echotemp.get()[tempprfid * PlusePoint + pluseid] = 
+							echotemp.get()[tempprfid * PlusePoint + pluseid]
+							+ std::complex<float>(h_echo[freqi].x, h_echo[freqi].y);
+						echoMasktemp(tempprfid, pluseid)=echoMasktemp(tempprfid, pluseid) + 1;
 					}
+
 					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;
 					}
@@ -915,7 +944,7 @@ ErrorCode RTPCProcessCls::RTPCMainProcess_GPU( )
 				}
 
 			}
-			
+			echoMaskImg.saveImage(echoMasktemp, startprfid, 0, 1);
 			this->EchoSimulationData->saveEchoArr(echotemp, startprfid, templine);
 
 		}
@@ -938,7 +967,7 @@ ErrorCode RTPCProcessCls::RTPCMainProcess_GPU( )
 	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_localangle);			FreeCUDADevice(d_localangle); //11
 
 	FreeCUDAHost(h_RstX);				FreeCUDADevice(d_RstX);
 	FreeCUDAHost(h_RstY);				FreeCUDADevice(d_RstY);
@@ -960,6 +989,65 @@ ErrorCode RTPCProcessCls::RTPCMainProcess_GPU( )
 	return ErrorCode::SUCCESS;
 }
 
+std::shared_ptr<SatelliteOribtNode[]> RTPCProcessCls::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 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)
 {

SimulationSAR/RTPCProcessCls.h

@@ -63,14 +63,19 @@ private:
 	QString TaskFileName;
 	QString tmpfolderPath;
 
+	QString OutEchoMaskPath;
 
 public:
 	ErrorCode Process(long num_thread); // 处理
 private: // 处理流程
 	ErrorCode InitParams();// 1. 初始化参数
 	ErrorCode DEMPreprocess(); // 2. 裁剪DEM范围
+	ErrorCode InitEchoMaskArray();
 	//ErrorCode RTPCMainProcess(long num_thread);
 	ErrorCode RTPCMainProcess_GPU();
+
+	std::shared_ptr<SatelliteOribtNode[]> getSatelliteOribtNodes(double prf_time, double dt, bool antflag, long double imageStarttime);
+
 private:
 	QString demxyzPath;
 	QString demmaskPath;

SimulationSAR/TBPImageAlgCls.cpp

@@ -42,7 +42,7 @@ void CreatePixelXYZ(std::shared_ptr<EchoL0Dataset> echoL0ds, QString outPixelXYZ
 		Eigen::MatrixXd demy = xyzRaster.getData(0, startcolidx, prfcount, tempechocol, 2);
 		Eigen::MatrixXd demz = xyzRaster.getData(0, startcolidx, prfcount, tempechocol, 3);
 
-
+		#pragma omp parallel for 
 		for (long i = 0; i < prfcount; i++) {
 
 			double Px = 0;
@@ -271,12 +271,17 @@ ErrorCode TBPImageAlgCls::ProcessGPU()
 				antpy.get()[anti] = Pys.get()[anti + startechoid];
 				antpz.get()[anti] = Pzs.get()[anti + startechoid];
 			}
-			TBPImageGPUAlg(antpx, antpy, antpz,
+			TBPImageGPUAlg(
+				antpx, antpy, antpz,
 				img_x, img_y, img_z,
 				echoArr, imgArr,
-				freq, fs, Rnear, Rfar,
+				freq, fs,
+				Rnear, Rfar,
 				tempimgBlockline, colCount, 
 				tempechoBlockline, PlusePoints );
+
+
+
 		}
 		this->L1ds->saveImageRaster(imgArr, startimgrowid, tempimgBlockline);
 	}
@@ -319,6 +324,12 @@ void TBPImageGPUAlg(std::shared_ptr<float> antPx, std::shared_ptr<float> antPy,
 	cuComplex* h_imgArr = (cuComplex*)mallocCUDAHost(sizeof(cuComplex) * rowcount * colcount);
 	cuComplex* d_imgArr = (cuComplex*)mallocCUDADevice( sizeof(cuComplex) * rowcount * colcount);
 
+
+	float* h_R=(float*)mallocCUDAHost(sizeof(float) * rowcount * colcount);
+	float* d_R=(float*)mallocCUDADevice(sizeof(float) * rowcount * colcount);
+
+
+
 	// ³õʼ»¯
 	// ÌìÏßλÖÃ
 	for (long i = 0; i < prfcount; i++) {
@@ -360,15 +371,18 @@ void TBPImageGPUAlg(std::shared_ptr<float> antPx, std::shared_ptr<float> antPy,
 	HostToDevice(h_imgx, d_imgx, sizeof(float) * rowcount * colcount);
 	HostToDevice(h_imgy, d_imgy, sizeof(float) * rowcount * colcount);
 	HostToDevice(h_imgz, d_imgz, sizeof(float) * rowcount * colcount);
+	HostToDevice(h_R, d_R, sizeof(float) * rowcount * colcount);
 
 	HostToDevice(h_echoArr, d_echoArr, sizeof(cuComplex) * prfcount * freqcount);
 	HostToDevice(h_imgArr, d_imgArr, sizeof(cuComplex) * rowcount * colcount);
- 
+
+
+
 
 	for (long prfid = 0; prfid < prfcount; prfid++) {
 		CUDATBPImage(
 			d_antPx,d_antPy,d_antPz,
-			d_imgx,d_imgy,d_imgz,
+			d_imgx,d_imgy,d_imgz, d_R,
 			d_echoArr,
 			d_imgArr,
 			freq, fs, Rnear, Rfar,
@@ -383,8 +397,8 @@ void TBPImageGPUAlg(std::shared_ptr<float> antPx, std::shared_ptr<float> antPy,
 
 	for (long i = 0; i < rowcount; i++) {
 		for (long j = 0; j < colcount; j++) {
-			imgArr.get()[i * colcount + j] = std::complex<double>(h_imgArr[i * colcount + j].x,
+			imgArr.get()[i * colcount + j] = imgArr.get()[i * colcount + j]+ std::complex<float>(h_imgArr[i * colcount + j].x,
-																  h_imgArr[i * colcount + j].y);
+																				h_imgArr[i * colcount + j].y);
 		}
 	}
 
@@ -392,11 +406,14 @@ void TBPImageGPUAlg(std::shared_ptr<float> antPx, std::shared_ptr<float> antPy,
 	FreeCUDAHost(h_antPx);		FreeCUDADevice(d_antPx);
 	FreeCUDAHost(h_antPy);		FreeCUDADevice(d_antPy);
 	FreeCUDAHost(h_antPz);		FreeCUDADevice(d_antPz);
+
 	FreeCUDAHost(h_imgx);		FreeCUDADevice(d_imgx);
 	FreeCUDAHost(h_imgy);		FreeCUDADevice(d_imgy);
 	FreeCUDAHost(h_imgz);		FreeCUDADevice(d_imgz);
+
 	FreeCUDAHost(h_echoArr);	FreeCUDADevice(d_echoArr);
 	FreeCUDAHost(h_imgArr);		FreeCUDADevice(d_imgArr);
+	FreeCUDAHost(h_R);			FreeCUDADevice(d_R);
 
 }