sarsimulation01.cpp

来自「C语言实现机载SAR成像」· C++ 代码 · 共 165 行

#include <iostream.h>
#include <fstream.h>
#include <math.h>
#include <stdlib.h>
#include <cstdlib>
#include <iomanip.h>
#define ARRAYNUM 2048
#define TARGETNUM 3
#define pi 3.14159
void main()
{
	//定义常数：光速,载频及波长
	double c=3e8;
	double fc=1e9;
	double lambda;
	lambda=c/fc;
	//定义目标区域
	double xmin=0,xmax=50;
	double yc=10000,y0=500;
	double v=100,H=5000;
	double R0=sqrt(yc*yc+H*H);
	//定义天线参数
	double D=4;
	double Lsar=lambda*R0/D;
	double Tsar=Lsar/v;
	//慢时间域参数
	double Ka=-2*v*v/lambda/R0;
	double Ba=Ka*Tsar;
	Ba=Ba>0?Ba:(-Ba);
	double PRF=Ba;
	double PRT=1/PRF;
	double ds=PRT;
	//慢时间域采样数
	int Nslow=int((xmax-xmin+Lsar)/v/ds);
	Nslow=Nslow+1;
	Nslow=int(log(Nslow)/log(2))+1;	
	int temp=1;
	for(int i=1;i<=Nslow;i++)
	{
		temp=temp*2;
	}
	Nslow=temp;
	//慢时间域离散采样序列
	double tp=((xmax+Lsar/2)/v-(xmin-Lsar/2)/v)/Nslow;
	double sn[ARRAYNUM];
	for(i=0;i<Nslow;i++)
	{
		sn[i]=(xmin-Lsar/2)/v+i*tp;
	}
	PRT=(xmax-xmin+Lsar)/v/Nslow;
	PRF=1/PRT;
	ds=PRT;
	//定义快时间域参数
	double Tr=5e-6;
	double Br=30e6;
	double Kr=Br/Tr;
	double Fsr=3*Br;
	double dt=1/Fsr;
	double Rmin=sqrt((yc-y0)*(yc-y0)+H*H);
	double Rmax=sqrt((yc+y0)*(yc+y0)+H*H+(Lsar/2)*(Lsar/2));
	//快时间域离散采样序列
	int Nfast=int(2*(Rmax-Rmin)/c/dt+Tr/dt)+1;
	Nfast=int(log(Nfast)/log(2))+1;	
	temp=1;
	for(i=1;i<=Nfast;i++)
	{
		temp=temp*2;
	}
	Nfast=temp;
	tp=(2*Rmax/c+Tr-2*Rmin/c)/Nfast;
	double tm[ARRAYNUM];
	for(i=0;i<Nfast;i++)
	{
		tm[i]=2*Rmin/c+i*tp;
	}
	dt=(2*Rmax/c+Tr-2*Rmin/c)/Nfast;
	Fsr=1/dt;
	//分辨率
	double DY=c/2/Br;
	double DX=D/2;
	//目标参量
	int Ntarget=3;
	double Ptarget[3][3]={{xmin,yc,1},{xmin,yc+10*DY,1},{xmin+10*DX,yc+20*DY,1}};
	int K=Ntarget;
	int N=Nslow;
	int M=Nfast;
	double T[TARGETNUM][3];
	for(i=0;i<3;i++)
		for(int j=0;j<3;j++)
		{
			T[i][j]=Ptarget[i][j];
		}
	//定义回波模型
	static double Srnmi[ARRAYNUM][ARRAYNUM],Srnmj[ARRAYNUM][ARRAYNUM],Dslow[ARRAYNUM],Dfast[ARRAYNUM][ARRAYNUM];
	for(i=0;i<N;i++)
		for(int j=0;j<M;j++)
		{
			Srnmi[i][j]=0;
			Srnmj[i][j]=0;
			Dfast[i][j]=0;
			Dslow[i]=0;
		}
    static double R[ARRAYNUM],tau[ARRAYNUM];
    static double Ptm[ARRAYNUM][ARRAYNUM],Ptau[ARRAYNUM][ARRAYNUM],PR[ARRAYNUM][ARRAYNUM],phase[ARRAYNUM][ARRAYNUM];
	for(int k=0;k<K;k++)
	{
		double sigma=T[k][2];
		for(i=0;i<Nslow;i++)
		{
			Dslow[i]=sn[i]*v-T[k][0];
			R[i]=sqrt(Dslow[i]*Dslow[i]+T[k][1]*T[k][1]+H*H);
			tau[i]=2*R[i]/c;
			Dslow[i]=Dslow[i]>0?Dslow[i]:(-Dslow[i]);
			if(Dslow[i]<Lsar/2)
				Dslow[i]=1;
			else
				Dslow[i]=0;
		}

		for(i=0;i<N;i++)
			for(int j=0;j<M;j++)
			{
				Ptm[i][j]=tm[j];
				Ptau[i][j]=tau[i];
				PR[i][j]=R[i];
				Dfast[i][j]=Ptm[i][j]-Ptau[i][j];
				phase[i][j]=pi*Kr*Dfast[i][j]*Dfast[i][j]-(4*pi/lambda)*PR[i][j];
				if((Dfast[i][j]>0)&&(Dfast[i][j])<Tr)
					Dfast[i][j]=1;
				else
					Dfast[i][j]=0;
			}
		for(i=0;i<N;i++)
			for(int j=0;j<M;j++)
			{
				Srnmi[i][j]=Srnmi[i][j]+sigma*cos(phase[i][j])*Dfast[i][j]*Dslow[i];
				Srnmj[i][j]=Srnmj[i][j]+sigma*sin(phase[i][j])*Dfast[i][j]*Dslow[i];
			}
	}
	ofstream out_stream;
	//将回波的实部存入sidata.dat中
	out_stream.open("sidata.dat");
	for(i=0;i<N;i++)
	{
		for(int j=0;j<M;j++)
		{
			out_stream<<setw(5)<<Srnmi[i][j]<<' ';
		}
		out_stream<<endl;
	}
	out_stream.close();
	ofstream out_stream1;
	//将回波的虚部部存入sjdata.dat中
	out_stream1.open("sjdata.dat");
	for(i=0;i<N;i++)
	{
		for(int j=0;j<M;j++)
		{
			out_stream1<<setw(5)<<Srnmj[i][j]<<' ';
		}
		out_stream1<<endl;
	}
	out_stream1.close();
	cout<<"End of output"<<endl;
}