M_map(五)

一、圓形區域的畫圖

1、

clear all

LATLIMS=[14 22];
LONLIMS=[108 118];%南海邊界范圍
m_proj('miller','lon',LONLIMS,'lat',LATLIMS);%投影模式
m_coast('patch',[0.1 0.7 0.4]);%綠色填充
m_grid('box','fancy','tickdir','in');%沒有網格，邊框相間，%m_line(lon,lat,'linewi',2,'color','r','linestyle',':');控制線條格式，點畫線還是直線

lon=112:1:116;
lat=16:1:20;
m_line(lon,lat,'linewi',2,'color','r');%線寬，2；顏色

[X,Y]=m_ll2xy(117,21);
line(X,Y,'marker','.','markersize',24','color','r')%畫點

2、距離

clear all
%distance用法
arclen=distance([37,0],[35,0])%返回兩點間的相對球心的角度，以度為單位
d=distance([37,0],[35,0],6371)% [緯度，經度] [緯度，經度] [半徑]
D=(arclen/180)*pi*6371
%m_map中函數
dist=m_lldist([0 0],[35 37])%[經度 經度] [緯度 緯度]

3、

%畫出矩形散點圖和圓形圖
clear all

LATLIMS=[14 22];
LONLIMS=[108 118];%南海邊界范圍
m_proj('lambert','lon',LONLIMS,'lat',LATLIMS);%投影模式
m_coast('patch',[0.1 0.7 0.4]);%綠色填充
m_grid('box','fancy','tickdir','in');%沒有網格，邊框相間，%m_line(lon,lat,'linewi',2,'color','r','linestyle',':');控制線條格式，點畫線還是直線

load EDH_south_sea_2008
load coordi_south_sea_2008 

m_range_ring(114.0623,17.9532,[1e2:1e2:3e2],'linewi',2,'color','b');%紅色300km范圍圓圈

% 矩形點陣
range_lat=4:24;%21N和15N對應的位置下標
range_lon=20:40;%111.5E和116.5E對應的下標
for i=1:length(range_lon)
    for j=1:length(range_lat)
        [X,Y]=m_ll2xy(lon_south_sea(range_lon(i)),lat_south_sea(range_lat(j)));%化為x,y坐標
        line(X,Y,'marker','.','markersize',10,'color','k')%畫點
        hold on
    end
end

%離散圓
[X0,Y0]=m_ll2xy(114.0623,17.9532);%圓心化為x,y坐標
line(X0,Y0,'marker','.','markersize',15,'color','r');%畫圓心
DIST=m_lldist([114 115],[18 18])%經度加1度，增加的距離
R=300;%300km
[X1,Y1]=m_ll2xy(114.0623-R/DIST,17.9532);%找到300km的一個點

r=sqrt((X0-X1)^2+(Y0-Y1)^2);%地圖距離到圖上距離轉換
r=linspace(0,r,10);%等分成9份
for i=1:10
    theta=0:pi/45:2*pi;%360度，91個點。0到2pi
    x=X0+r(i)*cos(theta); %（X0,Y0）圓心
    y=Y0+r(i)*sin(theta);
    plot(x,y,'.','markersize',6,'color','r')
end

　　

4、pcolor

clear all
n =18;
r = (0:n)'/n;
theta = pi*(-n:n)/n;
X = r*cos(theta);
Y = r*sin(theta);
C = r*cos(2*theta);
pcolor(X,Y,C)
axis equal tight
colorbar

figure
load PCOLOR %南海坐標和波導高度數據
colormap('jet');
shading flat;%平滑方式
gca=pcolor(Plg,Plt,EDH_south_sea)
set(gca, 'LineStyle','none');%去除網格
axis equal tight %按比例展示
colorbar %顏色條

clear all

LATLIMS=[2 22];
LONLIMS=[105 121];%南海邊界范圍
m_proj('miller','lon',LONLIMS,'lat',LATLIMS);%投影模式


load coordi_south_sea_2008 %載入南海經緯度信息
load EDH_south_sea_2008
for i=1:length(lon_south_sea) %經度
    [X,Y]=m_ll2xy(lon_south_sea(i),12);%化為x,y坐標,12是緯度的平均值
    lon(i)=X;
end
for i=1:length(lat_south_sea) %經度
    [X,Y]=m_ll2xy(113,lat_south_sea(i));%化為x,y坐標，113是經度的平均值
    lat(i)=Y;
end

[Plg,Plt]=meshgrid(lon,lat);%經緯度 形成網格

colormap('jet');
shading flat;%平滑方式
gca=pcolor(Plg,Plt,EDH_south_sea)
set(gca, 'LineStyle','none');%去除網格
axis equal tight %按比例展示
colorbar %顏色條

m_coast('patch',[0.1 0.7 0.4]);%綠色填充
m_grid('box','fancy','tickdir','in');　

左右兩幅圖對比，左邊是將經緯度換為xy值后畫的圖，更改繪圖模式，利用pcolor畫的，看出上邊邊緣有的吻合不是太好；右圖是直接利用m_pcolor畫的，吻合度更好

相當於載入了一幅地圖，先將其經緯度化為正常的坐標，然后再在上面畫圖。

5、

clear all
n =6;
r = (0:n)'/n;%0到6，半徑上均分的數
theta = pi*(-n:n)/n;%將整個圓分成了13分。
X = r*cos(theta);
Y = r*sin(theta);
C = r*cos(2*theta);
pcolor(X,Y,C)
axis equal tight
colorbar

角度逆時針轉

半徑分成了6份。

 6、

%溫度繪圖，只讀取南海數據繪圖
clear all
load EDH_south_sea_2008
load coordi_south_sea_2008 
LATLIMS=[14 22];
LONLIMS=[108 118];%南海邊界范圍
m_proj('lambert','lon',LONLIMS,'lat',LATLIMS);%投影模式

hold on %一定要有，否則地圖會被覆蓋

%% 
%將經緯度換成坐標軸
for i=1:length(lon_south_sea) %經度
    [X,Y]=m_ll2xy(lon_south_sea(i),18);%化為x,y坐標,12是緯度的平均值
    lon(i)=X;
end
for i=1:length(lat_south_sea) %經度
    [X,Y]=m_ll2xy(114,lat_south_sea(i));%化為x,y坐標，113是經度的平均值
    lat(i)=Y;
end

%% 
%圓形點圖數據
%m_range_ring(114.0623,17.9532,[3e2],'linewi',2,'color','b');%紅色300km范圍圓圈
%離散圓
[X0,Y0]=m_ll2xy(114.0623,17.9532);%圓心化為x,y坐標
%[X0,Y0]=m_ll2xy(114.0623,5);%圓心化為x,y坐標
%line(X0,Y0,'marker','.','markersize',15,'color','r');%畫圓心
DIST=m_lldist([114 115],[18 18]);%經度加1度，增加的距離
R=300;%300km
[X1,Y1]=m_ll2xy(114.0623-R/DIST,17.9532);%找到300km的一個點
%[X1,Y1]=m_ll2xy(114.0623-R/DIST,5);%找到300km的一個點
theta=0:pi/45:2*pi;%360度，91個點。0到2pi
r=sqrt((X0-X1)^2+(Y0-Y1)^2);%地圖距離到圖上距離轉換
r=linspace(0,r,10);%等分成10份
%theta=0:pi/45:2*pi;%360度，91個點。0到2pi

%% 
%定義裝位置和數值的空矩陣
circle_x=zeros(length(r),length(theta));
circle_y=zeros(length(r),length(theta));
circle_EDH=zeros(length(r),length(theta));
%% 
%N*2D方式，一共生成了91條線，每條線上有10個數據，
for i=1:length(r)
    x=X0+r(i)*cos(theta); %（X0,Y0）圓心，得到想要獲取的EDH的位置（x,y）
    y=Y0+r(i)*sin(theta);
    for j=1:length(theta)
        [min_lon,lon_index]=min(abs(lon-x(j)));%得到最匹配的x軸的值
        [min_lat,lat_index]=min(abs(lat-y(j)));%得到最匹配的y軸的值
        circle_EDH(i,j)=EDH_south_sea(lat_index,lon_index);%存儲蒸發波導高度
    end
    circle_x(i,:)=x;
    circle_y(i,:)=y;
    %plot(x,y,'.','markersize',6,'color','r')
end
%% 
%畫圖
colormap('jet');
shading flat;%平滑方式
gca=pcolor(circle_x,circle_y,circle_EDH)
set(gca, 'LineStyle','none');%去除網格
axis equal tight %按比例展示
% colorbar %顏色條

c=colorbar('eastoutside','fontsize',12)
%caxis([0,14])%顏色條范圍限制
title( '2008年南海平均蒸發波導高度','fontsize',15);%[a b c d]確定繪圖區域的位置, [a b]為繪圖區域左下點的坐標。
m_coast('patch',[0.1 0.7 0.4]);%綠色填充
m_grid('box','fancy','tickdir','in');         %c，d分別為繪圖區域的寬和高。

這樣可以看到，靠近岸的地方溫度低，靠近海的地方溫度高。

7、獲取路徑損失數據圖

%獲取圓形區域路徑損失數據
clear all;
% load testUsrDef.mat
tic
freq = 8000;    % in MHz頻率
thetabw = 2;    % in deg 3分貝波束寬度
thetae = 0; % in deg 仰角度數
polrz = 1;  % 1 for horizontal polarization 水平極化, 2 for vertical polarization 極化方式
tx_height = 3;  % tx_antenna height, in m，天線高度
range = 300; % calcu range, in km，水平傳播距離
zmax_user = 100;    % maximum height (max desired calculation height), in m，電波限高
rx_height=3;%接收天線高度
delx = 100; % range step, in m，距離步進值迭代設置
delz = 0.1;    % altitude step, in m，高度步進值
%%
%%地形設置
edge_range = [ 20 30 50 ];   % array of edge range, in km，地形限制
edge_height = [ 5 20 7 ];    % array of edge heights at edge ranges, in m，地形高度；？？？無效？
terrain_type = 1;   % tarrain type, 1 for no terrain case, 2 for terrain case 
interp_type = 2;    % edge interp type, 2 for linear, 3 for cubic spline
backward = 1;   % 1 for 1-way SSPE, 2 for 2-way SSPE

%% 海洋表面參數
ground_type = 1;    % ground type, 1 for PEC, 2 for mixed
epsilon = 69.13;    % dielectric constant 介電常數
sigma = 7.146;  % conductivity 電導率
%% 
%預留數據空間
load CIRCLE_EDH  % 10*91 91次循環，91條線，10個蒸發波導高度
circle_path_loss=zeros(10*range+1,length(circle_EDH));
circle_pl_x=zeros(10*range+1,length(circle_EDH));
circle_pl_y=zeros(10*range+1,length(circle_EDH));

%% 波導條件設置
% duct_type = 5;  % duct type, 1 for standard atmosphere, 3 for surface-based duct, 5 for evaporation duct
% duct_M =  [360 0];%必須帶0；
% duct_height = [15 0]; %必須帶0
% duct_range = 1; %？？？水平不均勻性的設置

% duct_type = 5;  % duct type, 1 for standard atmosphere, 3 for surface-based duct, 5 for evaporation duct
% duct_M =  [360 0 ;300 0;344 0];
% duct_height = [15 0;18 0; 20 0]; 
% duct_range =[0 100 200]; %水平不均勻性的設置,0km處蒸發波導高度15m,100km處蒸發波導高度5m,200km處蒸發波導高度20m.
for i=1:length(circle_EDH)
    duct_type = 5;  % duct type, 1 for standard atmosphere, 3 for surface-based duct, 5 for evaporation duct
    
    duct_M=zeros(10,2);
    
    duct_height = zeros(10,2);
    duct_height(:,1)=circle_EDH(:,i);%獲取第i列（線）波導高度
    duct_range =linspace(0,300,10); %水平不均勻性的設置,0km處蒸發波導高度15m,100km處蒸發波導高度5m,200km處蒸發波導高度20m.
    %計算路徑損失
    [umat, path_loss, prop_fact, free_space_loss, range_vec, z_user, z, stopflag] =  SSPE_function(freq, ...
        thetabw, thetae, polrz, tx_height, range, zmax_user, edge_range, edge_height, ...
        duct_type, duct_M, duct_height, duct_range, terrain_type, interp_type, backward, ...
        ground_type, epsilon, sigma, delx, delz);
    circle_path_loss(:,i)=path_loss(rx_height*10,:);
    
end
%計算坐標
for i=1:length(circle_EDH)
   circle_pl_x(:,i)=linspace(circle_x(1,i),circle_x(10,i),10*range+1);
   circle_pl_y(:,i)=linspace(circle_y(1,i),circle_y(10,i),10*range+1);
end
toc
save CIRCLE_PATH_LOSS circle_path_loss circle_pl_x circle_pl_y

8、圓形路徑損失圖

%圓形路徑損失繪制
clear all
LATLIMS=[14 22];
LONLIMS=[108 118];%南海邊界范圍
m_proj('lambert','lon',LONLIMS,'lat',LATLIMS);%投影模式
hold on %一定要有，否則地圖會被覆蓋
%% 
load CIRCLE_PATH_LOSS
%% 
%畫圖
colormap('jet');
shading flat;%平滑方式
gca=pcolor(circle_pl_x,circle_pl_y,circle_path_loss);
set(gca, 'LineStyle','none');%去除網格
axis equal tight ;%按比例展示
%colorbar %顏色條

c=colorbar('eastoutside','fontsize',12);
caxis([100,220]);%顏色條范圍限制
m_coast('patch',[0.1 0.7 0.4]);%綠色填充
m_grid('box','fancy','tickdir','in');         %c，d分別為繪圖區域的寬和高。