KJ135
/
RasterProcessTool
mirror of http://172.16.0.12:5000/LAMPSARToolSoftware/RasterProcessTool.git
4
0
Fork 0
Code Issues Packages Projects Releases Wiki Activity

【测试】废弃成像的单精度版本，因为单精度计算精度存在问题

pull/8/head
chenzenghui 2025-03-06 20:20:12 +08:00
parent 29d0eae76d
commit eb9fe9e1c1
3 changed files with 396 additions and 659 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

263
Toolbox/SimulationSARTool/SimulationSAR/BPBasic0_CUDA.cu
View File

@ -185,10 +185,7 @@ void bpBasic0CUDA(GPUDATA& data, int flag,double* h_R) {
data.im_final data.im_final
//,d_R //,d_R
); );
PrintLasterError("processPulseKernel");
if (ii % 1000==0) {
printfinfo("\rPRF(%f %) %d / %d\t\t\t\t",(ii*100.0/data.Np), ii,data.Np);
}
} }
//FreeCUDADevice(d_R); //FreeCUDADevice(d_R);
@ -280,264 +277,6 @@ void BPBasic0(GPUDATA& h_data)
/** 单精度代码**********************************************************************************************/
__global__ void phaseCompensationKernel_single(cufftComplex* phdata, const float* Freq, float r, int K, int Na) {
int freqIdx = blockIdx.x * blockDim.x + threadIdx.x;
int pulseIdx = blockIdx.y * blockDim.y + threadIdx.y;
if (freqIdx >= K || pulseIdx >= Na) return;
int idx = pulseIdx * K + freqIdx;
float phase = 4 * PI * Freq[freqIdx] * r / c;
float cos_phase = cosf(phase);
float sin_phase = sinf(phase);
cufftComplex ph = phdata[idx];
float new_real = ph.x * cos_phase - ph.y * sin_phase;
float new_imag = ph.x * sin_phase + ph.y * cos_phase;
phdata[idx] = make_cuComplex(new_real, new_imag);
}
__global__ void fftshiftKernel_single(cufftComplex* data, int Nfft, int Np) {
int pulse = blockIdx.x * blockDim.x + threadIdx.x;
if (pulse >= Np) return;
int half = Nfft / 2;
for (int i = 0; i < half; ++i) {
cufftComplex temp = data[pulse * Nfft + i];
data[pulse * Nfft + i] = data[pulse * Nfft + i + half];
data[pulse * Nfft + i + half] = temp;
}
}
__global__ void processPulseKernel_single(
long prfid,
int nx, int ny,
const float* x_mat, const float* y_mat, const float* z_mat,
float AntX, float AntY, float AntZ,
float R0, float minF,
const cufftComplex* rc_pulse,
const float r_start, const float dr, const int nR,
cufftComplex* im_final
) {
//
long long idx = blockIdx.x * blockDim.x + threadIdx.x;
long long pixelcount = nx * ny;
if (idx >= pixelcount) return;
//printf("processPulseKernel start!!\n");
//if (x >= nx || y >= ny) return;
//int idx = x * ny + y;
float dx = AntX - x_mat[idx];
float dy = AntY - y_mat[idx];
float dz = AntZ - z_mat[idx];
//printf("processPulseKernel xmat !!\n");
float R = sqrtf(dx * dx + dy * dy + dz * dz);
float dR = R - R0;
if (dR < r_start || dR >= (r_start + dr * (nR - 1))) return;
// Linear interpolation
float pos = (dR - r_start) / dr;
int index = (int)floorf(pos);
float weight = pos - index;
if (index < 0 || index >= nR - 1) return;
cufftComplex rc_low = rc_pulse[prfid * nR + index];
cufftComplex rc_high = rc_pulse[prfid * nR + index + 1];
cufftComplex rc_interp;
rc_interp.x = rc_low.x * (1 - weight) + rc_high.x * weight;
rc_interp.y = rc_low.y * (1 - weight) + rc_high.y * weight;
// Phase correction
float phase = 4 * PI * minF / c * dR;
float cos_phase = cosf(phase);
float sin_phase = sinf(phase);
cufftComplex phCorr;
phCorr.x = rc_interp.x * cos_phase - rc_interp.y * sin_phase;
phCorr.y = rc_interp.x * sin_phase + rc_interp.y * cos_phase;
// Accumulate
im_final[idx].x += phCorr.x;
im_final[idx].y += phCorr.y;
//printf("r_start=%e;dr=%e;nR=%d\n", r_start, dr, nR);
}
void bpBasic0CUDA_single(GPUDATA_single& data, int flag, float* h_R)
{
// Phase compensation
if (flag == 1) {
dim3 block(16, 16);
dim3 grid((data.K + 15) / 16, (data.Np + 15) / 16);
phaseCompensationKernel_single << <grid, block >> > (data.phdata, data.Freq, data.R0, data.K, data.Np);
PrintLasterError("bpBasic0CUDA Phase compensation");
//data.R0 = data.r; // 假设data.r已正确设置
}
// FFT处理
cufftHandle plan;
cufftPlan1d(&plan, data.Nfft, CUFFT_C2C, data.Np);
cufftExecC2C(plan, data.phdata, data.phdata, CUFFT_INVERSE);
cufftDestroy(plan);
// FFT移位
dim3 blockShift(256);
dim3 gridShift((data.Np + 255) / 256);
fftshiftKernel_single << <gridShift, blockShift >> > (data.phdata, data.Nfft, data.Np);
PrintLasterError("bpBasic0CUDA Phase FFT Process");
printfinfo("fft finished!!\n");
// 图像重建
float r_start = data.r_vec[0];
float dr = (data.r_vec[data.Nfft - 1] - r_start) / (data.Nfft - 1);
printfinfo("dr = %f\n", dr);
long pixelcount = data.nx * data.ny;
long grid_size = (pixelcount + BLOCK_SIZE - 1) / BLOCK_SIZE;
printfinfo("grid finished!!\n");
//float* d_R = (float*)mallocCUDADevice(sizeof(float) * data.nx * data.ny);
printfinfo("r_start=%e;dr=%e;nR=%d\n", r_start, dr, data.Nfft);
printfinfo("BPimage .....\n");
for (long ii = 0; ii < data.Np; ++ii) {
processPulseKernel_single << <grid_size, BLOCK_SIZE >> > (
ii,
data.nx, data.ny,
data.x_mat, data.y_mat, data.z_mat,
data.AntX[ii], data.AntY[ii], data.AntZ[ii],
data.R0, data.minF[ii],
data.phdata,
r_start, dr, data.Nfft,
data.im_final
//,d_R
);
PrintLasterError("processPulseKernel");
if (ii % 1000 == 0) {
printfinfo("\rPRF(%f %) %d / %d\t\t\t\t", (ii * 100.0 / data.Np), ii, data.Np);
}
}
//FreeCUDADevice(d_R);
PrintLasterError("bpBasic0CUDA Phase BPimage Process finished!!");
}
void initGPUData_single(GPUDATA_single& h_data, GPUDATA_single& d_data)
{
d_data.AntX = h_data.AntX; //(double*)mallocCUDADevice(sizeof(double) * h_data.Np);
d_data.AntY = h_data.AntY;//(double*)mallocCUDADevice(sizeof(double) * h_data.Np);
d_data.AntZ = h_data.AntZ;// (double*)mallocCUDADevice(sizeof(double) * h_data.Np);
d_data.minF = h_data.minF;// (double*)mallocCUDADevice(sizeof(double) * h_data.Np);
d_data.x_mat = (float*)mallocCUDADevice(sizeof(float) * h_data.nx * h_data.ny);
d_data.y_mat = (float*)mallocCUDADevice(sizeof(float) * h_data.nx * h_data.ny);
d_data.z_mat = (float*)mallocCUDADevice(sizeof(float) * h_data.nx * h_data.ny);
d_data.r_vec = h_data.r_vec;// (double*)mallocCUDADevice(sizeof(double) * h_data.Nfft);
d_data.Freq = (float*)mallocCUDADevice(sizeof(float) * h_data.Nfft);
d_data.phdata = (cuComplex*)mallocCUDADevice(sizeof(cuComplex) * h_data.Nfft * h_data.Np);
d_data.im_final = (cuComplex*)mallocCUDADevice(sizeof(cuComplex) * h_data.nx * h_data.ny);
//HostToDevice(h_data.AntX, d_data.AntX,sizeof(double) * h_data.Np);
//HostToDevice(h_data.AntY, d_data.AntY,sizeof(double) * h_data.Np);
//HostToDevice(h_data.AntZ, d_data.AntZ,sizeof(double) * h_data.Np);
//HostToDevice(h_data.minF, d_data.minF,sizeof(double) * h_data.Np);
HostToDevice(h_data.x_mat, d_data.x_mat, sizeof(float) * h_data.nx * h_data.ny); printf("image X Copy finished!!!\n");
HostToDevice(h_data.y_mat, d_data.y_mat, sizeof(float) * h_data.nx * h_data.ny); printf("image Y Copy finished!!!\n");
HostToDevice(h_data.z_mat, d_data.z_mat, sizeof(float) * h_data.nx * h_data.ny); printf("image Z Copy finished!!!\n");
HostToDevice(h_data.Freq, d_data.Freq, sizeof(float) * h_data.Nfft);
//HostToDevice(h_data.r_vec, d_data.r_vec, sizeof(double) * h_data.Nfft);
HostToDevice(h_data.phdata, d_data.phdata, sizeof(cuComplex) * h_data.Nfft * h_data.Np); printf("image echo Copy finished!!!\n");
HostToDevice(h_data.im_final, d_data.im_final, sizeof(cuComplex) * h_data.nx * h_data.ny); printf("image data Copy finished!!!\n");
// 拷贝标量参数
d_data.Nfft = h_data.Nfft;
d_data.K = h_data.K;
d_data.Np = h_data.Np;
d_data.nx = h_data.nx;
d_data.ny = h_data.ny;
d_data.R0 = h_data.R0;
d_data.deltaF = h_data.deltaF;
}
void freeGPUData_single(GPUDATA_single& d_data)
{
//FreeCUDADevice((d_data.AntX));
//FreeCUDADevice((d_data.AntY));
//FreeCUDADevice((d_data.AntZ));
//FreeCUDADevice((d_data.minF));
FreeCUDADevice((d_data.x_mat));
FreeCUDADevice((d_data.y_mat));
FreeCUDADevice((d_data.z_mat));
//FreeCUDADevice((d_data.r_vec));
FreeCUDADevice((d_data.Freq));
FreeCUDADevice((d_data.phdata));
FreeCUDADevice((d_data.im_final));
}
void freeHostData_single(GPUDATA_single& h_data)
{
//FreeCUDAHost((h_data.AntX));
//FreeCUDAHost((h_data.AntY));
//FreeCUDAHost((h_data.AntZ));
FreeCUDAHost((h_data.minF));
//FreeCUDAHost((h_data.x_mat));
//FreeCUDAHost((h_data.y_mat));
//FreeCUDAHost((h_data.z_mat));
FreeCUDAHost((h_data.r_vec));
FreeCUDAHost((h_data.Freq));
FreeCUDAHost((h_data.phdata));
FreeCUDAHost((h_data.im_final));
}
void BPBasic0_single(GPUDATA_single& h_data)
{
GPUDATA_single d_data;
initGPUData_single(h_data, d_data);
bpBasic0CUDA_single(d_data, 0);
DeviceToHost(h_data.im_final, d_data.im_final, sizeof(cuComplex) * h_data.nx * h_data.ny);
freeGPUData_single(d_data);
}
//int main() { //int main() {
// GPUDATA h_data, d_data; // GPUDATA h_data, d_data;
// //

11
Toolbox/SimulationSARTool/SimulationSAR/BPBasic0_CUDA.cuh
View File

@ -55,11 +55,12 @@ extern "C" {
extern "C" { extern "C" {
void bpBasic0CUDA_single(GPUDATA_single& data, int flag, float* h_R = nullptr); // 函数废弃，因为这个成像精度不够
void initGPUData_single(GPUDATA_single& h_data, GPUDATA_single& d_data); //void bpBasic0CUDA_single(GPUDATA_single& data, int flag, float* h_R = nullptr);
void freeGPUData_single(GPUDATA_single& d_data); //void initGPUData_single(GPUDATA_single& h_data, GPUDATA_single& d_data);
void freeHostData_single(GPUDATA_single& d_data); //void freeGPUData_single(GPUDATA_single& d_data);
void BPBasic0_single(GPUDATA_single& h_data); //void freeHostData_single(GPUDATA_single& d_data);
//void BPBasic0_single(GPUDATA_single& h_data);
}; };

777
Toolbox/SimulationSARTool/UnitTestMain.cpp
View File

@ -404,393 +404,393 @@ void testBpImage() {
} }
//
//
void testSimualtionEchoPoint_single() { //void testSimualtionEchoPoint_single() {
GPUDATA_single h_data{}; // GPUDATA_single h_data{};
/** 1. 轨道 **************************************************************************************************/ // /** 1. 轨道 **************************************************************************************************/
qDebug() << u8"1.轨道文件读取中。。。"; // qDebug() << u8"1.轨道文件读取中。。。";
QString inGPSPath = u8"C:\\Users\\30453\\Desktop\\script\\data\\GF3_Simulation.gpspos.data"; // QString inGPSPath = u8"C:\\Users\\30453\\Desktop\\script\\data\\GF3_Simulation.gpspos.data";
long gpspoints = gdalImage(inGPSPath).height; // long gpspoints = gdalImage(inGPSPath).height;
std::shared_ptr<SatelliteAntPos> antpos = SatelliteAntPosOperator::readAntPosFile(inGPSPath, gpspoints); // std::shared_ptr<SatelliteAntPos> antpos = SatelliteAntPosOperator::readAntPosFile(inGPSPath, gpspoints);
h_data.AntX = (float*)mallocCUDAHost(sizeof(float) * gpspoints); // h_data.AntX = (float*)mallocCUDAHost(sizeof(float) * gpspoints);
h_data.AntY = (float*)mallocCUDAHost(sizeof(float) * gpspoints); // h_data.AntY = (float*)mallocCUDAHost(sizeof(float) * gpspoints);
h_data.AntZ = (float*)mallocCUDAHost(sizeof(float) * gpspoints); // h_data.AntZ = (float*)mallocCUDAHost(sizeof(float) * gpspoints);
//
for (long i = 0; i < gpspoints; i = i + 1) { // for (long i = 0; i < gpspoints; i = i + 1) {
h_data.AntX[i] = antpos.get()[i].Px; // h_data.AntX[i] = antpos.get()[i].Px;
h_data.AntY[i] = antpos.get()[i].Py; // h_data.AntY[i] = antpos.get()[i].Py;
h_data.AntZ[i] = antpos.get()[i].Pz; // h_data.AntZ[i] = antpos.get()[i].Pz;
} // }
//gpspoints = gpspoints / 2; // //gpspoints = gpspoints / 2;
qDebug() << "GPS points Count:\t" << gpspoints; // qDebug() << "GPS points Count:\t" << gpspoints;
qDebug() << "PRF 0:\t " << QString("%1,%2,%3").arg(h_data.AntX[0]).arg(h_data.AntY[0]).arg(h_data.AntZ[0]); // qDebug() << "PRF 0:\t " << QString("%1,%2,%3").arg(h_data.AntX[0]).arg(h_data.AntY[0]).arg(h_data.AntZ[0]);
qDebug() << "PRF " << gpspoints - 1 << ":\t " << QString("%1,%2,%3") // qDebug() << "PRF " << gpspoints - 1 << ":\t " << QString("%1,%2,%3")
.arg(h_data.AntX[gpspoints - 1]) // .arg(h_data.AntX[gpspoints - 1])
.arg(h_data.AntY[gpspoints - 1]) // .arg(h_data.AntY[gpspoints - 1])
.arg(h_data.AntZ[gpspoints - 1]); // .arg(h_data.AntZ[gpspoints - 1]);
/** 2. 成像网格 **************************************************************************************************/ // /** 2. 成像网格 **************************************************************************************************/
qDebug() << u8"轨道文件读取结束\n2.成像网格读取。。。"; // qDebug() << u8"轨道文件读取结束\n2.成像网格读取。。。";
QString demxyzPath = u8"C:\\Users\\30453\\Desktop\\script\\已修改GF3场景\\data\\demxyz.bin"; // QString demxyzPath = u8"C:\\Users\\30453\\Desktop\\script\\已修改GF3场景\\data\\demxyz.bin";
gdalImage demgridimg(demxyzPath); // gdalImage demgridimg(demxyzPath);
long dem_rowCount = demgridimg.height; // long dem_rowCount = demgridimg.height;
long dem_ColCount = demgridimg.width; // long dem_ColCount = demgridimg.width;
std::shared_ptr<float> demX = readDataArr<float>(demgridimg, 0, 0, dem_rowCount, dem_ColCount, 1, GDALREADARRCOPYMETHOD::VARIABLEMETHOD); // std::shared_ptr<float> demX = readDataArr<float>(demgridimg, 0, 0, dem_rowCount, dem_ColCount, 1, GDALREADARRCOPYMETHOD::VARIABLEMETHOD);
std::shared_ptr<float> demY = readDataArr<float>(demgridimg, 0, 0, dem_rowCount, dem_ColCount, 2, GDALREADARRCOPYMETHOD::VARIABLEMETHOD); // std::shared_ptr<float> demY = readDataArr<float>(demgridimg, 0, 0, dem_rowCount, dem_ColCount, 2, GDALREADARRCOPYMETHOD::VARIABLEMETHOD);
std::shared_ptr<float> demZ = readDataArr<float>(demgridimg, 0, 0, dem_rowCount, dem_ColCount, 3, GDALREADARRCOPYMETHOD::VARIABLEMETHOD); // std::shared_ptr<float> demZ = readDataArr<float>(demgridimg, 0, 0, dem_rowCount, dem_ColCount, 3, GDALREADARRCOPYMETHOD::VARIABLEMETHOD);
//
h_data.x_mat = (float*)mallocCUDAHost(sizeof(float) * dem_rowCount * dem_ColCount); // h_data.x_mat = (float*)mallocCUDAHost(sizeof(float) * dem_rowCount * dem_ColCount);
h_data.y_mat = (float*)mallocCUDAHost(sizeof(float) * dem_rowCount * dem_ColCount); // h_data.y_mat = (float*)mallocCUDAHost(sizeof(float) * dem_rowCount * dem_ColCount);
h_data.z_mat = (float*)mallocCUDAHost(sizeof(float) * dem_rowCount * dem_ColCount); // h_data.z_mat = (float*)mallocCUDAHost(sizeof(float) * dem_rowCount * dem_ColCount);
//
for (long i = 0; i < dem_rowCount; i++) { // for (long i = 0; i < dem_rowCount; i++) {
for (long j = 0; j < dem_ColCount; j++) { // for (long j = 0; j < dem_ColCount; j++) {
h_data.x_mat[i * dem_ColCount + j] = demX.get()[i * dem_ColCount + j]; // h_data.x_mat[i * dem_ColCount + j] = demX.get()[i * dem_ColCount + j];
h_data.y_mat[i * dem_ColCount + j] = demY.get()[i * dem_ColCount + j]; // h_data.y_mat[i * dem_ColCount + j] = demY.get()[i * dem_ColCount + j];
h_data.z_mat[i * dem_ColCount + j] = demZ.get()[i * dem_ColCount + j]; // h_data.z_mat[i * dem_ColCount + j] = demZ.get()[i * dem_ColCount + j];
} // }
} // }
qDebug() << "dem row Count:\t" << dem_rowCount; // qDebug() << "dem row Count:\t" << dem_rowCount;
qDebug() << "dem col Count:\t" << dem_ColCount; // qDebug() << "dem col Count:\t" << dem_ColCount;
//
qDebug() << u8"成像网格读取结束"; // qDebug() << u8"成像网格读取结束";
/** 3. 频率网格 **************************************************************************************************/ // /** 3. 频率网格 **************************************************************************************************/
qDebug() << u8"3.处理频率"; // qDebug() << u8"3.处理频率";
float centerFreq = 5.3e9; // float centerFreq = 5.3e9;
float bandwidth = 40e6; // float bandwidth = 40e6;
long freqpoints = 2048; // long freqpoints = 2048;
std::shared_ptr<float> freqlist(getFreqPoints_mallocHost_single(centerFreq - bandwidth / 2, centerFreq + bandwidth / 2, freqpoints), FreeCUDAHost); // std::shared_ptr<float> freqlist(getFreqPoints_mallocHost_single(centerFreq - bandwidth / 2, centerFreq + bandwidth / 2, freqpoints), FreeCUDAHost);
std::shared_ptr<float> minF(new float[gpspoints], delArrPtr); // std::shared_ptr<float> minF(new float[gpspoints], delArrPtr);
for (long i = 0; i < gpspoints; i++) { // for (long i = 0; i < gpspoints; i++) {
minF.get()[i] = freqlist.get()[0]; // minF.get()[i] = freqlist.get()[0];
} // }
h_data.deltaF = bandwidth / (freqpoints - 1); // h_data.deltaF = bandwidth / (freqpoints - 1);
qDebug() << "start Freq:\t" << centerFreq - bandwidth / 2; // qDebug() << "start Freq:\t" << centerFreq - bandwidth / 2;
qDebug() << "end Freq:\t" << centerFreq + bandwidth / 2; // qDebug() << "end Freq:\t" << centerFreq + bandwidth / 2;
qDebug() << "freq points:\t" << freqpoints; // qDebug() << "freq points:\t" << freqpoints;
qDebug() << "delta freq:\t" << freqlist.get()[1] - freqlist.get()[0]; // qDebug() << "delta freq:\t" << freqlist.get()[1] - freqlist.get()[0];
qDebug() << u8"频率结束"; // qDebug() << u8"频率结束";
/** 4. 初始化回波 **************************************************************************************************/ // /** 4. 初始化回波 **************************************************************************************************/
qDebug() << u8"4.初始化回波"; // qDebug() << u8"4.初始化回波";
std::shared_ptr<cuComplex> phdata(createEchoPhase_mallocHost(gpspoints, freqpoints), FreeCUDAHost); // std::shared_ptr<cuComplex> phdata(createEchoPhase_mallocHost(gpspoints, freqpoints), FreeCUDAHost);
qDebug() << "Azimuth Points:\t" << gpspoints; // qDebug() << "Azimuth Points:\t" << gpspoints;
qDebug() << "Range Points:\t" << freqpoints; // qDebug() << "Range Points:\t" << freqpoints;
qDebug() << u8"初始化回波结束"; // qDebug() << u8"初始化回波结束";
/** 5. 初始化图像 **************************************************************************************************/ // /** 5. 初始化图像 **************************************************************************************************/
qDebug() << u8"5.初始化图像"; // qDebug() << u8"5.初始化图像";
h_data.im_final = (cuComplex*)mallocCUDAHost(sizeof(cuComplex) * dem_rowCount * dem_ColCount); // h_data.im_final = (cuComplex*)mallocCUDAHost(sizeof(cuComplex) * dem_rowCount * dem_ColCount);
qDebug() << "Azimuth Points:\t" << gpspoints; // qDebug() << "Azimuth Points:\t" << gpspoints;
qDebug() << "Range Points:\t" << freqpoints; // qDebug() << "Range Points:\t" << freqpoints;
qDebug() << u8"初始化图像结束"; // qDebug() << u8"初始化图像结束";
//
/** 6. 模型参数初始化 **************************************************************************************************/ // /** 6. 模型参数初始化 **************************************************************************************************/
qDebug() << u8"6.模型参数初始化"; // qDebug() << u8"6.模型参数初始化";
//
h_data.nx = dem_ColCount; // h_data.nx = dem_ColCount;
h_data.ny = dem_rowCount; // h_data.ny = dem_rowCount;
//
h_data.Np = gpspoints; // h_data.Np = gpspoints;
h_data.Freq = freqlist.get(); // h_data.Freq = freqlist.get();
h_data.minF = minF.get(); // h_data.minF = minF.get();
h_data.Nfft = freqpoints; // h_data.Nfft = freqpoints;
h_data.K = h_data.Nfft; // h_data.K = h_data.Nfft;
h_data.phdata = phdata.get(); // h_data.phdata = phdata.get();
h_data.R0 = 900000; // h_data.R0 = 900000;
qDebug() << u8"模型参数结束"; // qDebug() << u8"模型参数结束";
//
//
/** 7. 目标 **************************************************************************************************/ // /** 7. 目标 **************************************************************************************************/
float Tx = -2029086.618142, Ty = 4139594.934504, Tz = 4392846.782027; // float Tx = -2029086.618142, Ty = 4139594.934504, Tz = 4392846.782027;
float Tslx = -2029086.618142, Tsly = 4139594.934504, Tslz = 4392846.782027; // float Tslx = -2029086.618142, Tsly = 4139594.934504, Tslz = 4392846.782027;
float p1 = 1, p2 = 0, p3 = 0, p4 = 0, p5 = 0, p6 = 0; // float p1 = 1, p2 = 0, p3 = 0, p4 = 0, p5 = 0, p6 = 0;
//
/** 7. 构建回波 **************************************************************************************************/ // /** 7. 构建回波 **************************************************************************************************/
GPUDATA_single d_data; // GPUDATA_single d_data;
initGPUData_single(h_data, d_data); // initGPUData_single(h_data, d_data);
//
RFPCProcess_single(Tx, Ty, Tz, // RFPCProcess_single(Tx, Ty, Tz,
Tslx, Tsly, Tslz, // 目标的坡面向量 // Tslx, Tsly, Tslz, // 目标的坡面向量
p1, p2, p3, p4, p5, p6, // p1, p2, p3, p4, p5, p6,
d_data); // d_data);
//
//
/** 8. 展示回波 **************************************************************************************************/ // /** 8. 展示回波 **************************************************************************************************/
{ // {
DeviceToHost(h_data.phdata, d_data.phdata, sizeof(cuComplex) * d_data.Np * d_data.Nfft); // DeviceToHost(h_data.phdata, d_data.phdata, sizeof(cuComplex) * d_data.Np * d_data.Nfft);
cuComplex* h_ifftphdata = (cuComplex*)mallocCUDAHost(sizeof(cuComplex) * d_data.Np * d_data.Nfft); // cuComplex* h_ifftphdata = (cuComplex*)mallocCUDAHost(sizeof(cuComplex) * d_data.Np * d_data.Nfft);
cuComplex* d_ifftphdata = (cuComplex*)mallocCUDADevice(sizeof(cuComplex) * d_data.Np * d_data.Nfft); // cuComplex* d_ifftphdata = (cuComplex*)mallocCUDADevice(sizeof(cuComplex) * d_data.Np * d_data.Nfft);
CUDAIFFT(d_data.phdata, d_ifftphdata, d_data.Np, d_data.Nfft, d_data.Nfft); // CUDAIFFT(d_data.phdata, d_ifftphdata, d_data.Np, d_data.Nfft, d_data.Nfft);
FFTShift1D(d_ifftphdata, d_data.Np, d_data.Nfft); // FFTShift1D(d_ifftphdata, d_data.Np, d_data.Nfft);
DeviceToHost(h_ifftphdata, d_ifftphdata, sizeof(cuComplex) * d_data.Np * d_data.Nfft); // DeviceToHost(h_ifftphdata, d_ifftphdata, sizeof(cuComplex) * d_data.Np * d_data.Nfft);
std::shared_ptr<std::complex<double>> ifftdata(new std::complex<double>[d_data.Np * d_data.Nfft], delArrPtr); // std::shared_ptr<std::complex<double>> ifftdata(new std::complex<double>[d_data.Np * d_data.Nfft], delArrPtr);
{ // {
for (long i = 0; i < d_data.Np; i++) { // for (long i = 0; i < d_data.Np; i++) {
for (long j = 0; j < d_data.Nfft; j++) { // for (long j = 0; j < d_data.Nfft; j++) {
ifftdata.get()[i * d_data.Nfft + j] = // ifftdata.get()[i * d_data.Nfft + j] =
std::complex<double>( // std::complex<double>(
h_ifftphdata[i * d_data.Nfft + j].x, // h_ifftphdata[i * d_data.Nfft + j].x,
h_ifftphdata[i * d_data.Nfft + j].y); // h_ifftphdata[i * d_data.Nfft + j].y);
} // }
} // }
} // }
testOutComplexDoubleArr(QString("echo_ifft_single.bin"), ifftdata.get(), d_data.Np, d_data.Nfft); // testOutComplexDoubleArr(QString("echo_ifft_single.bin"), ifftdata.get(), d_data.Np, d_data.Nfft);
//
FreeCUDADevice(d_ifftphdata); // FreeCUDADevice(d_ifftphdata);
FreeCUDAHost(h_ifftphdata); // FreeCUDAHost(h_ifftphdata);
//
} // }
/** 9. 成像 **************************************************************************************************/ // /** 9. 成像 **************************************************************************************************/
//
// 计算maxWr（需要先计算deltaF） // // 计算maxWr（需要先计算deltaF）
float deltaF = h_data.deltaF; // 从输入参数获取 // float deltaF = h_data.deltaF; // 从输入参数获取
float maxWr = 299792458.0f / (2.0f * deltaF); // float maxWr = 299792458.0f / (2.0f * deltaF);
qDebug() << "maxWr :\t" << maxWr; // qDebug() << "maxWr :\t" << maxWr;
float* r_vec_host = (float*)mallocCUDAHost(sizeof(float) * h_data.Nfft);// new float[data.Nfft]; // float* r_vec_host = (float*)mallocCUDAHost(sizeof(float) * h_data.Nfft);// new float[data.Nfft];
const float step = maxWr / h_data.Nfft; // const float step = maxWr / h_data.Nfft;
const float start = -1 * h_data.Nfft / 2.0f * step; // const float start = -1 * h_data.Nfft / 2.0f * step;
printf("nfft=%d\n", h_data.Nfft); // printf("nfft=%d\n", h_data.Nfft);
//
for (int i = 0; i < h_data.Nfft; ++i) { // for (int i = 0; i < h_data.Nfft; ++i) {
r_vec_host[i] = start + i * step; // r_vec_host[i] = start + i * step;
} // }
//
h_data.r_vec = r_vec_host; // h_data.r_vec = r_vec_host;
d_data.r_vec = h_data.r_vec; // d_data.r_vec = h_data.r_vec;
//
qDebug() << "r_vec [0]:\t" << h_data.r_vec[0]; // qDebug() << "r_vec [0]:\t" << h_data.r_vec[0];
qDebug() << "r_vec step:\t" << step; // qDebug() << "r_vec step:\t" << step;
//
bpBasic0CUDA_single(d_data, 0); // bpBasic0CUDA_single(d_data, 0);
//
DeviceToHost(h_data.im_final, d_data.im_final, sizeof(cuComplex) * d_data.nx * d_data.ny); // DeviceToHost(h_data.im_final, d_data.im_final, sizeof(cuComplex) * d_data.nx * d_data.ny);
//
{ // {
DeviceToHost(h_data.phdata, d_data.phdata, sizeof(cuComplex) * d_data.Np * d_data.Nfft); // DeviceToHost(h_data.phdata, d_data.phdata, sizeof(cuComplex) * d_data.Np * d_data.Nfft);
std::shared_ptr<std::complex<double>> ifftdata(new std::complex<double>[d_data.Np * d_data.Nfft], delArrPtr); // std::shared_ptr<std::complex<double>> ifftdata(new std::complex<double>[d_data.Np * d_data.Nfft], delArrPtr);
{ // {
for (long i = 0; i < d_data.Np; i++) { // for (long i = 0; i < d_data.Np; i++) {
for (long j = 0; j < d_data.Nfft; j++) { // for (long j = 0; j < d_data.Nfft; j++) {
ifftdata.get()[i * d_data.Nfft + j] = // ifftdata.get()[i * d_data.Nfft + j] =
std::complex<double>( // std::complex<double>(
h_data.phdata[i * d_data.Nfft + j].x, // h_data.phdata[i * d_data.Nfft + j].x,
h_data.phdata[i * d_data.Nfft + j].y); // h_data.phdata[i * d_data.Nfft + j].y);
} // }
} // }
} // }
testOutComplexDoubleArr(QString("echo_ifft_BPBasic_single.bin"), ifftdata.get(), d_data.Np, d_data.Nfft); // testOutComplexDoubleArr(QString("echo_ifft_BPBasic_single.bin"), ifftdata.get(), d_data.Np, d_data.Nfft);
//
//
//
std::shared_ptr<std::complex<double>> im_finals(new std::complex<double>[d_data.nx * d_data.ny], delArrPtr); // std::shared_ptr<std::complex<double>> im_finals(new std::complex<double>[d_data.nx * d_data.ny], delArrPtr);
{ // {
for (long i = 0; i < d_data.ny; i++) { // for (long i = 0; i < d_data.ny; i++) {
for (long j = 0; j < d_data.nx; j++) { // for (long j = 0; j < d_data.nx; j++) {
im_finals.get()[i * d_data.nx + j] = std::complex<double>( // im_finals.get()[i * d_data.nx + j] = std::complex<double>(
h_data.im_final[i * d_data.nx + j].x, // h_data.im_final[i * d_data.nx + j].x,
h_data.im_final[i * d_data.nx + j].y); // h_data.im_final[i * d_data.nx + j].y);
} // }
} // }
} // }
testOutComplexDoubleArr(QString("im_finals_single.bin"), im_finals.get(), d_data.ny, d_data.nx); // testOutComplexDoubleArr(QString("im_finals_single.bin"), im_finals.get(), d_data.ny, d_data.nx);
} // }
} //}
//
//
void testBpImage_single() { //void testBpImage_single() {
//
GPUDATA_single h_data{}; // GPUDATA_single h_data{};
/** 0. 轨道 **************************************************************************************************/ // /** 0. 轨道 **************************************************************************************************/
QString echoPath = "D:\\Programme\\vs2022\\RasterMergeTest\\LAMPCAE_SCANE-all-scane\\GF3_Simulation.xml"; // QString echoPath = "D:\\Programme\\vs2022\\RasterMergeTest\\LAMPCAE_SCANE-all-scane\\GF3_Simulation.xml";
std::shared_ptr<EchoL0Dataset> echoL0ds(new EchoL0Dataset); // std::shared_ptr<EchoL0Dataset> echoL0ds(new EchoL0Dataset);
echoL0ds->Open(echoPath); // echoL0ds->Open(echoPath);
qDebug() << u8"加载回拨文件：\t" << echoPath; // qDebug() << u8"加载回拨文件：\t" << echoPath;
//
/** 1. 轨道 **************************************************************************************************/ // /** 1. 轨道 **************************************************************************************************/
qDebug() << u8"1.轨道文件读取中。。。"; // qDebug() << u8"1.轨道文件读取中。。。";
QString inGPSPath = echoL0ds->getGPSPointFilePath(); // QString inGPSPath = echoL0ds->getGPSPointFilePath();
long gpspoints = gdalImage(inGPSPath).height; // long gpspoints = gdalImage(inGPSPath).height;
std::shared_ptr<SatelliteAntPos> antpos = SatelliteAntPosOperator::readAntPosFile(inGPSPath, gpspoints); // std::shared_ptr<SatelliteAntPos> antpos = SatelliteAntPosOperator::readAntPosFile(inGPSPath, gpspoints);
h_data.AntX = (float*)mallocCUDAHost(sizeof(float) * gpspoints); // h_data.AntX = (float*)mallocCUDAHost(sizeof(float) * gpspoints);
h_data.AntY = (float*)mallocCUDAHost(sizeof(float) * gpspoints); // h_data.AntY = (float*)mallocCUDAHost(sizeof(float) * gpspoints);
h_data.AntZ = (float*)mallocCUDAHost(sizeof(float) * gpspoints); // h_data.AntZ = (float*)mallocCUDAHost(sizeof(float) * gpspoints);
//
for (long i = 0; i < gpspoints; i = i + 1) { // for (long i = 0; i < gpspoints; i = i + 1) {
h_data.AntX[i] = antpos.get()[i].Px; // h_data.AntX[i] = antpos.get()[i].Px;
h_data.AntY[i] = antpos.get()[i].Py; // h_data.AntY[i] = antpos.get()[i].Py;
h_data.AntZ[i] = antpos.get()[i].Pz; // h_data.AntZ[i] = antpos.get()[i].Pz;
} // }
//gpspoints = gpspoints / 2; // //gpspoints = gpspoints / 2;
qDebug() << "GPS points Count:\t" << gpspoints; // qDebug() << "GPS points Count:\t" << gpspoints;
qDebug() << "PRF 0:\t " << QString("%1,%2,%3").arg(h_data.AntX[0]).arg(h_data.AntY[0]).arg(h_data.AntZ[0]); // qDebug() << "PRF 0:\t " << QString("%1,%2,%3").arg(h_data.AntX[0]).arg(h_data.AntY[0]).arg(h_data.AntZ[0]);
qDebug() << "PRF " << gpspoints - 1 << ":\t " << QString("%1,%2,%3") // qDebug() << "PRF " << gpspoints - 1 << ":\t " << QString("%1,%2,%3")
.arg(h_data.AntX[gpspoints - 1]) // .arg(h_data.AntX[gpspoints - 1])
.arg(h_data.AntY[gpspoints - 1]) // .arg(h_data.AntY[gpspoints - 1])
.arg(h_data.AntZ[gpspoints - 1]); // .arg(h_data.AntZ[gpspoints - 1]);
/** 2. 成像网格 **************************************************************************************************/ // /** 2. 成像网格 **************************************************************************************************/
qDebug() << u8"轨道文件读取结束\n2.成像网格读取。。。"; // qDebug() << u8"轨道文件读取结束\n2.成像网格读取。。。";
QString demxyzPath = u8"D:\\Programme\\vs2022\\RasterMergeTest\\simulationData\\demdataset\\demxyz.bin"; // QString demxyzPath = u8"D:\\Programme\\vs2022\\RasterMergeTest\\simulationData\\demdataset\\demxyz.bin";
gdalImage demgridimg(demxyzPath); // gdalImage demgridimg(demxyzPath);
long dem_rowCount = demgridimg.height; // long dem_rowCount = demgridimg.height;
long dem_ColCount = demgridimg.width; // long dem_ColCount = demgridimg.width;
std::shared_ptr<double> demX = readDataArr<double>(demgridimg, 0, 0, dem_rowCount, dem_ColCount, 1, GDALREADARRCOPYMETHOD::VARIABLEMETHOD); // std::shared_ptr<double> demX = readDataArr<double>(demgridimg, 0, 0, dem_rowCount, dem_ColCount, 1, GDALREADARRCOPYMETHOD::VARIABLEMETHOD);
std::shared_ptr<double> demY = readDataArr<double>(demgridimg, 0, 0, dem_rowCount, dem_ColCount, 2, GDALREADARRCOPYMETHOD::VARIABLEMETHOD); // std::shared_ptr<double> demY = readDataArr<double>(demgridimg, 0, 0, dem_rowCount, dem_ColCount, 2, GDALREADARRCOPYMETHOD::VARIABLEMETHOD);
std::shared_ptr<double> demZ = readDataArr<double>(demgridimg, 0, 0, dem_rowCount, dem_ColCount, 3, GDALREADARRCOPYMETHOD::VARIABLEMETHOD); // std::shared_ptr<double> demZ = readDataArr<double>(demgridimg, 0, 0, dem_rowCount, dem_ColCount, 3, GDALREADARRCOPYMETHOD::VARIABLEMETHOD);
//
h_data.x_mat = (float*)mallocCUDAHost(sizeof(float) * dem_rowCount * dem_ColCount); // h_data.x_mat = (float*)mallocCUDAHost(sizeof(float) * dem_rowCount * dem_ColCount);
h_data.y_mat = (float*)mallocCUDAHost(sizeof(float) * dem_rowCount * dem_ColCount); // h_data.y_mat = (float*)mallocCUDAHost(sizeof(float) * dem_rowCount * dem_ColCount);
h_data.z_mat = (float*)mallocCUDAHost(sizeof(float) * dem_rowCount * dem_ColCount); // h_data.z_mat = (float*)mallocCUDAHost(sizeof(float) * dem_rowCount * dem_ColCount);
//
for (long i = 0; i < dem_rowCount; i++) { // for (long i = 0; i < dem_rowCount; i++) {
for (long j = 0; j < dem_ColCount; j++) { // for (long j = 0; j < dem_ColCount; j++) {
h_data.x_mat[i * dem_ColCount + j] = demX.get()[i * dem_ColCount + j]; // h_data.x_mat[i * dem_ColCount + j] = demX.get()[i * dem_ColCount + j];
h_data.y_mat[i * dem_ColCount + j] = demY.get()[i * dem_ColCount + j]; // h_data.y_mat[i * dem_ColCount + j] = demY.get()[i * dem_ColCount + j];
h_data.z_mat[i * dem_ColCount + j] = demZ.get()[i * dem_ColCount + j]; // h_data.z_mat[i * dem_ColCount + j] = demZ.get()[i * dem_ColCount + j];
} // }
} // }
qDebug() << "dem row Count:\t" << dem_rowCount; // qDebug() << "dem row Count:\t" << dem_rowCount;
qDebug() << "dem col Count:\t" << dem_ColCount; // qDebug() << "dem col Count:\t" << dem_ColCount;
//
qDebug() << u8"成像网格读取结束"; // qDebug() << u8"成像网格读取结束";
/** 3. 频率网格 **************************************************************************************************/ // /** 3. 频率网格 **************************************************************************************************/
qDebug() << u8"3.处理频率"; // qDebug() << u8"3.处理频率";
float centerFreq = echoL0ds->getCenterFreq(); // float centerFreq = echoL0ds->getCenterFreq();
float bandwidth = echoL0ds->getBandwidth(); // float bandwidth = echoL0ds->getBandwidth();
size_t freqpoints = echoL0ds->getPlusePoints(); // size_t freqpoints = echoL0ds->getPlusePoints();
std::shared_ptr<float> freqlist(getFreqPoints_mallocHost_single(centerFreq - bandwidth / 2, centerFreq + bandwidth / 2, freqpoints), FreeCUDAHost); // std::shared_ptr<float> freqlist(getFreqPoints_mallocHost_single(centerFreq - bandwidth / 2, centerFreq + bandwidth / 2, freqpoints), FreeCUDAHost);
std::shared_ptr<float> minF(new float[gpspoints], delArrPtr); // std::shared_ptr<float> minF(new float[gpspoints], delArrPtr);
for (long i = 0; i < gpspoints; i++) { // for (long i = 0; i < gpspoints; i++) {
minF.get()[i] = freqlist.get()[0]; // minF.get()[i] = freqlist.get()[0];
} // }
h_data.deltaF = bandwidth / (freqpoints - 1); // h_data.deltaF = bandwidth / (freqpoints - 1);
qDebug() << "start Freq:\t" << centerFreq - bandwidth / 2; // qDebug() << "start Freq:\t" << centerFreq - bandwidth / 2;
qDebug() << "end Freq:\t" << centerFreq + bandwidth / 2; // qDebug() << "end Freq:\t" << centerFreq + bandwidth / 2;
qDebug() << "freq points:\t" << freqpoints; // qDebug() << "freq points:\t" << freqpoints;
qDebug() << "delta freq:\t" << freqlist.get()[1] - freqlist.get()[0]; // qDebug() << "delta freq:\t" << freqlist.get()[1] - freqlist.get()[0];
qDebug() << u8"频率结束"; // qDebug() << u8"频率结束";
/** 4. 初始化回波 **************************************************************************************************/ // /** 4. 初始化回波 **************************************************************************************************/
qDebug() << u8"4.初始化回波"; // qDebug() << u8"4.初始化回波";
std::shared_ptr<std::complex<double>> echoData = echoL0ds->getEchoArr(); // std::shared_ptr<std::complex<double>> echoData = echoL0ds->getEchoArr();
size_t echosize = sizeof(cuComplex) * gpspoints * freqpoints; // size_t echosize = sizeof(cuComplex) * gpspoints * freqpoints;
qDebug() << "echo data size (byte) :\t" << echosize; // qDebug() << "echo data size (byte) :\t" << echosize;
h_data.phdata = (cuComplex*)mallocCUDAHost(echosize); // h_data.phdata = (cuComplex*)mallocCUDAHost(echosize);
shared_complexPtrToHostCuComplex(echoData.get(), h_data.phdata, gpspoints * freqpoints); // shared_complexPtrToHostCuComplex(echoData.get(), h_data.phdata, gpspoints * freqpoints);
//
qDebug() << "Azimuth Points:\t" << gpspoints; // qDebug() << "Azimuth Points:\t" << gpspoints;
qDebug() << "Range Points:\t" << freqpoints; // qDebug() << "Range Points:\t" << freqpoints;
qDebug() << u8"初始化回波结束"; // qDebug() << u8"初始化回波结束";
/** 5. 初始化图像 **************************************************************************************************/ // /** 5. 初始化图像 **************************************************************************************************/
qDebug() << u8"5.初始化图像"; // qDebug() << u8"5.初始化图像";
h_data.im_final = (cuComplex*)mallocCUDAHost(sizeof(cuComplex) * dem_rowCount * dem_ColCount); // h_data.im_final = (cuComplex*)mallocCUDAHost(sizeof(cuComplex) * dem_rowCount * dem_ColCount);
qDebug() << "Azimuth Points:\t" << gpspoints; // qDebug() << "Azimuth Points:\t" << gpspoints;
qDebug() << "Range Points:\t" << freqpoints; // qDebug() << "Range Points:\t" << freqpoints;
qDebug() << u8"初始化图像结束"; // qDebug() << u8"初始化图像结束";
//
/** 6. 模型参数初始化 **************************************************************************************************/ // /** 6. 模型参数初始化 **************************************************************************************************/
qDebug() << u8"6.模型参数初始化"; // qDebug() << u8"6.模型参数初始化";
//
h_data.nx = dem_ColCount; // h_data.nx = dem_ColCount;
h_data.ny = dem_rowCount; // h_data.ny = dem_rowCount;
//
h_data.Np = gpspoints; // h_data.Np = gpspoints;
h_data.Freq = freqlist.get(); // h_data.Freq = freqlist.get();
h_data.minF = minF.get(); // h_data.minF = minF.get();
h_data.Nfft = freqpoints; // h_data.Nfft = freqpoints;
h_data.K = h_data.Nfft; // h_data.K = h_data.Nfft;
//
h_data.R0 = echoL0ds->getRefPhaseRange(); // h_data.R0 = echoL0ds->getRefPhaseRange();
qDebug() << u8"模型参数结束"; // qDebug() << u8"模型参数结束";
//
//
/** 7. 目标 **************************************************************************************************/ // /** 7. 目标 **************************************************************************************************/
double Tx = -2029086.618142, Ty = 4139594.934504, Tz = 4392846.782027; // double Tx = -2029086.618142, Ty = 4139594.934504, Tz = 4392846.782027;
double Tslx = -2029086.618142, Tsly = 4139594.934504, Tslz = 4392846.782027; // double Tslx = -2029086.618142, Tsly = 4139594.934504, Tslz = 4392846.782027;
double p1 = 33.3, p2 = 0, p3 = 0, p4 = 0, p5 = 0, p6 = 0; // double p1 = 33.3, p2 = 0, p3 = 0, p4 = 0, p5 = 0, p6 = 0;
//
/** 7. 构建回波 **************************************************************************************************/ // /** 7. 构建回波 **************************************************************************************************/
GPUDATA_single d_data; // GPUDATA_single d_data;
initGPUData_single(h_data, d_data); // initGPUData_single(h_data, d_data);
//
/** 8. 展示回波 **************************************************************************************************/ // /** 8. 展示回波 **************************************************************************************************/
{ // {
DeviceToHost(h_data.phdata, d_data.phdata, sizeof(cuComplex) * d_data.Np * d_data.Nfft); // DeviceToHost(h_data.phdata, d_data.phdata, sizeof(cuComplex) * d_data.Np * d_data.Nfft);
cuComplex* h_ifftphdata = (cuComplex*)mallocCUDAHost(sizeof(cuComplex) * d_data.Np * d_data.Nfft); // cuComplex* h_ifftphdata = (cuComplex*)mallocCUDAHost(sizeof(cuComplex) * d_data.Np * d_data.Nfft);
cuComplex* d_ifftphdata = (cuComplex*)mallocCUDADevice(sizeof(cuComplex) * d_data.Np * d_data.Nfft); // cuComplex* d_ifftphdata = (cuComplex*)mallocCUDADevice(sizeof(cuComplex) * d_data.Np * d_data.Nfft);
CUDAIFFT(d_data.phdata, d_ifftphdata, d_data.Np, d_data.Nfft, d_data.Nfft); // CUDAIFFT(d_data.phdata, d_ifftphdata, d_data.Np, d_data.Nfft, d_data.Nfft);
FFTShift1D(d_ifftphdata, d_data.Np, d_data.Nfft); // FFTShift1D(d_ifftphdata, d_data.Np, d_data.Nfft);
DeviceToHost(h_ifftphdata, d_ifftphdata, sizeof(cuComplex) * d_data.Np * d_data.Nfft); // DeviceToHost(h_ifftphdata, d_ifftphdata, sizeof(cuComplex) * d_data.Np * d_data.Nfft);
std::shared_ptr<std::complex<double>> ifftdata(new std::complex<double>[d_data.Np * d_data.Nfft], delArrPtr); // std::shared_ptr<std::complex<double>> ifftdata(new std::complex<double>[d_data.Np * d_data.Nfft], delArrPtr);
{ // {
for (long i = 0; i < d_data.Np; i++) { // for (long i = 0; i < d_data.Np; i++) {
for (long j = 0; j < d_data.Nfft; j++) { // for (long j = 0; j < d_data.Nfft; j++) {
ifftdata.get()[i * d_data.Nfft + j] = // ifftdata.get()[i * d_data.Nfft + j] =
std::complex<double>( // std::complex<double>(
h_ifftphdata[i * d_data.Nfft + j].x, // h_ifftphdata[i * d_data.Nfft + j].x,
h_ifftphdata[i * d_data.Nfft + j].y); // h_ifftphdata[i * d_data.Nfft + j].y);
} // }
} // }
} // }
testOutComplexDoubleArr(QString("echo_ifft.bin"), ifftdata.get(), d_data.Np, d_data.Nfft); // testOutComplexDoubleArr(QString("echo_ifft_single.bin"), ifftdata.get(), d_data.Np, d_data.Nfft);
//
FreeCUDADevice(d_ifftphdata); // FreeCUDADevice(d_ifftphdata);
FreeCUDAHost(h_ifftphdata); // FreeCUDAHost(h_ifftphdata);
//
} // }
/** 9. 成像 **************************************************************************************************/ // /** 9. 成像 **************************************************************************************************/
//
// 计算maxWr（需要先计算deltaF） // // 计算maxWr（需要先计算deltaF）
double deltaF = h_data.deltaF; // 从输入参数获取 // double deltaF = h_data.deltaF; // 从输入参数获取
double maxWr = 299792458.0f / (2.0f * deltaF); // double maxWr = 299792458.0f / (2.0f * deltaF);
qDebug() << "maxWr :\t" << maxWr; // qDebug() << "maxWr :\t" << maxWr;
float* r_vec_host = (float*)mallocCUDAHost(sizeof(float) * h_data.Nfft);// new double[data.Nfft]; // float* r_vec_host = (float*)mallocCUDAHost(sizeof(float) * h_data.Nfft);// new double[data.Nfft];
const double step = maxWr / h_data.Nfft; // const double step = maxWr / h_data.Nfft;
const double start = -1 * h_data.Nfft / 2.0f * step; // const double start = -1 * h_data.Nfft / 2.0f * step;
printf("nfft=%d\n", h_data.Nfft); // printf("nfft=%d\n", h_data.Nfft);
//
for (int i = 0; i < h_data.Nfft; ++i) { // for (int i = 0; i < h_data.Nfft; ++i) {
r_vec_host[i] = start + i * step; // r_vec_host[i] = start + i * step;
} // }
//
h_data.r_vec = r_vec_host; // h_data.r_vec = r_vec_host;
d_data.r_vec = h_data.r_vec; // d_data.r_vec = h_data.r_vec;
//
qDebug() << "r_vec [0]:\t" << h_data.r_vec[0]; // qDebug() << "r_vec [0]:\t" << h_data.r_vec[0];
qDebug() << "r_vec step:\t" << step; // qDebug() << "r_vec step:\t" << step;
//
bpBasic0CUDA_single(d_data, 0); // bpBasic0CUDA_single(d_data, 0);
//
DeviceToHost(h_data.im_final, d_data.im_final, sizeof(cuComplex) * d_data.nx * d_data.ny); // DeviceToHost(h_data.im_final, d_data.im_final, sizeof(cuComplex) * d_data.nx * d_data.ny);
{ // {
DeviceToHost(h_data.phdata, d_data.phdata, sizeof(cuComplex) * d_data.Np * d_data.Nfft); // DeviceToHost(h_data.phdata, d_data.phdata, sizeof(cuComplex) * d_data.Np * d_data.Nfft);
std::shared_ptr<std::complex<double>> ifftdata(new std::complex<double>[d_data.Np * d_data.Nfft], delArrPtr); // std::shared_ptr<std::complex<double>> ifftdata(new std::complex<double>[d_data.Np * d_data.Nfft], delArrPtr);
{ // {
for (long i = 0; i < d_data.Np; i++) { // for (long i = 0; i < d_data.Np; i++) {
for (long j = 0; j < d_data.Nfft; j++) { // for (long j = 0; j < d_data.Nfft; j++) {
ifftdata.get()[i * d_data.Nfft + j] = // ifftdata.get()[i * d_data.Nfft + j] =
std::complex<double>( // std::complex<double>(
h_data.phdata[i * d_data.Nfft + j].x, // h_data.phdata[i * d_data.Nfft + j].x,
h_data.phdata[i * d_data.Nfft + j].y); // h_data.phdata[i * d_data.Nfft + j].y);
} // }
} // }
} // }
testOutComplexDoubleArr(QString("echo_ifft_BPBasic.bin"), ifftdata.get(), d_data.Np, d_data.Nfft); // testOutComplexDoubleArr(QString("echo_ifft_BPBasic_single.bin"), ifftdata.get(), d_data.Np, d_data.Nfft);
//
std::shared_ptr<std::complex<double>> im_finals(new std::complex<double>[d_data.nx * d_data.ny], delArrPtr); // std::shared_ptr<std::complex<double>> im_finals(new std::complex<double>[d_data.nx * d_data.ny], delArrPtr);
{ // {
for (long i = 0; i < d_data.ny; i++) { // for (long i = 0; i < d_data.ny; i++) {
for (long j = 0; j < d_data.nx; j++) { // for (long j = 0; j < d_data.nx; j++) {
im_finals.get()[i * d_data.nx + j] = std::complex<double>( // im_finals.get()[i * d_data.nx + j] = std::complex<double>(
h_data.im_final[i * d_data.nx + j].x, // h_data.im_final[i * d_data.nx + j].x,
h_data.im_final[i * d_data.nx + j].y); // h_data.im_final[i * d_data.nx + j].y);
} // }
} // }
} // }
testOutComplexDoubleArr(QString("im_finals.bin"), im_finals.get(), d_data.ny, d_data.nx); // testOutComplexDoubleArr(QString("im_finals_single.bin"), im_finals.get(), d_data.ny, d_data.nx);
} // }
//
//
//
} //}
//
@ -799,9 +799,6 @@ int main(int argc, char* argv[]) {
QApplication a(argc, argv); QApplication a(argc, argv);
//testSimualtionEchoPoint();
//testSimualtionEchoPoint_single();
void testBpImage_single();
return 0; return 0;
} }