- 积分
- 19
- 贡献
-
- 精华
- 在线时间
- 小时
- 注册时间
- 2017-4-16
- 最后登录
- 1970-1-1
|
登录后查看更多精彩内容~
您需要 登录 才可以下载或查看,没有帐号?立即注册
x
各位大神,可以帮忙看下我的程序哪出问题了吗?
我是要做雷达山体阻挡分析的,可每次运行程序都出不来结果,我的循环结构还可以改进吗?
%DEM read
DEM=imread('E:\MATLAB\DEM\Guangdong 90\DEM_Guangdong90m.tif');
lonDEM=109.999583817611++0.000833333332+(0:1:12002)*0.00083333333;
latDEM=25.0004171267299-0.00083333333/2-(0:1:6001)*0.00083333333;
%Radar site informationfor guangzhou
RDST.Name='guangzhou';
RDST.lon=113.3550;
RDST.lat=23.00389;
RDST.El=180.6;
h0=RDST.El/1000;%雷达高程
%雷达山体阻挡分析
Stzd=zeros(9,360,230);
Re=8500.0;
%elv=FxAgl;
elv=[0.5,1.5,2.4,3.4,4.3,6.0, 9.9,14.6,19.5];
Els=elv*pi/180;
r=(1:1:230)-0.5;
az=0:1:359;
for i=1:9
iEl=Els(i); % Do the first elvation scan analysis
for j=1:360 %interate azimuth angle
iAz=az(j);
for k=1:230; %interate the radar gate along the radar beam
kR=r(k);
h=sqrt(kR*kR+Re*Re+2*kR*Re*sin(iEl))-Re;
kH=h+h0; % Now get the elevation for kth radar gate
Rs=Re*asin(kR*cos(iEl)/(Re+kH)); %obtain the surface distance from the radar site and the point of kth radar gate
[iLon,iLat]=fnt_RangeAzToLonLat(RDST.lon,RDST.lat, iAz,Rs); %obtain the DEM at location (iLon,iLat)
iDEM=fnt_GetDEM(iLon,iLat);
%detLat=abs(Lat-ilat);
%detLon=abs(Lon-ilon);
%m=find(detLat==min(detLat));
%n=find(detLon==min(detLon));
%iDEM=DEM(m,n);
if kH< iDEM %block
Stzd(i,j,k:230)=1;
else % kH>iDEM not block
Stzd(i,j,k)=0;
end
end
end
end
|
|