1. M_MAP を使ってベクトルなどを描画

描画例

冬季の海氷密接度（shaded）、海氷ベクトル、海面気圧、海上風を 1枚に重ねた例。

スクリプト例

旧ページの内容をそのまま使いやすい形で残しています。コメントの赤字部分は、 何をしている処理かを説明する注釈です。

clear ; clc ; clf ; close all
%% make freezing season mean field ice, Uice, SLP and UVwnd

%------- read latlon info. ---------------　AMSRデータの緯度経度情報を読み込み
fid = fopen('/home/kmizobat/work/AMSR/AMSR-PS/IC0/latlon_low_NP','r');
latlon = fread(fid,[900,1800],'short');
lat = latlon(:,1:900)/100. ; lon = latlon(:,901:1800)/100.;

%% ------- read IBCAO flat file 10801x751 2-byte integer [65-90N, 0-360E]
ibf = ['/home/kmizobat/work/MAP/IBCAO/IBCAO_ver2_23_GEO_netCDF_2min.flat'];　予め用意したIBCAOのデータを読み込み
fidib=fopen(ibf,'r') ; ibd=fread(fidib,[10801,751],'short'); fclose(fidib);
ibd2=-ibd(5401:8101,:); % Extract
sprintf('%s','finish reading IBCAO')

iblon=zeros(2701,751) ; iblat=iblon;
for j=1:751
    iblat(:,j)=90-(j-1)*((90-65)/751);
    for i=1:2701
    iblon(i,:)=-180+(i-1)*(360/10801);
    end
end
for jk=1:751
    for ik=1:2701
        if(iblon(ik,jk) < 180 & iblon(ik,jk) >= 0)
            ibd2(ik,jk)=NaN;
        end
        if(iblat(ik,jk) > 81 | iblat(ik,jk) < 65)
            ibd2(ik,jk)=NaN;
        end
    end
end
sprintf('%s','finish modifying IBCAO')

%% loop start　ここからForループの始まり
for yr=2002:2009
  yr2=num2str(yr);
  icef=['/home/kmizobat/work/AMSR-kamo/winter/sim_raw/mat/FALLMEAN/' yr2 '.fallmean.flat'] %海氷密接度ファイル指定
  uicef=['/home/kmizobat/work/AMSR-kamo/winter/sim_raw/mat/FALLMEAN/' yr2 '.fallmean.mve.mat'] %氷速ベクトル東西成分のファイル指定
  vicef=['/home/kmizobat/work/AMSR-kamo/winter/sim_raw/mat/FALLMEAN/' yr2 '.fallmean.mvn.mat'] %氷速ベクトル南北成分のファイル指定
  uwndf=['/home/kmizobat/work/AMSR-kamo/winter/sim_raw/mat/FALLMEAN/uwnd.' yr2 '.fallmean.flat'] %風速ベクトル東西成分のファイル指定
  vwndf=['/home/kmizobat/work/AMSR-kamo/winter/sim_raw/mat/FALLMEAN/vwnd.' yr2 '.fallmean.flat'] %風速ベクトル南北成分のファイル指定
  slpf=['/home/kmizobat/work/AMSR-kamo/winter/sim_raw/mat/FALLMEAN/slp.' yr2 '.fallmean.flat'] %海面気圧のファイル指定
  ofl=['/home/kmizobat/work/AMSR-kamo/winter/sim_raw/mat/FALLMEAN/',yr2,'.frzup.jpg']; ofl %出力ファイル名の指定

  fid = fopen(icef,'r'); ice = fread(fid,[900,900],'float'); clear fid;
  load(uicef) ; load(vicef);
  fid = fopen(uwndf,'r'); uwnd = fread(fid,[192,94],'float'); clear fid;
  fid = fopen(vwndf,'r'); vwnd = fread(fid,[192,94],'float'); clear fid;
  fid = fopen(slpf,'r'); slp = fread(fid,[144,72],'float'); clear fid;

  %% modifying ice conc. data % 欠損データの除去
  for jk=1:900 ; for ik=1:900
      if(lon(ik,jk) < 180 & lon(ik,jk) > 0)
          ice(ik,jk)=NaN;
      end
      if(lat(ik,jk) > 81)
          ice(ik,jk)=NaN;
      end
  end ; end

  %% modifying icevec. data
  for jk=1:151 ; for ik=1:159
      if(mve(ik,jk) < -8000 )
          mve(ik,jk)=NaN;
      end
      if(mvn(ik,jk) < -8000 )
          mvn(ik,jk)=NaN;
      end
  end ; end

  %% modifying slp data
  plon=zeros(144,72) ; plat=zeros(144,72);
  for jj=1:72 ; for ii=1:144
      plon(ii,jj)=0+(ii-1)*2.5;
      plat(ii,jj)=-90+(jj-1)*2.5;
      if(plon(ii,jj) < 180 & plon(ii,jj) > 0)
          slp(ii,jj)=NaN;
      end
      if(plon(ii,jj) > -120 & plon(ii,jj) <= 0)
          slp(ii,jj)=NaN;
      end
  end ; end
  ind=find(plon >= 180) ; plon(ind)=plon(ind)-360; clear ind;

  %% modifying uvwnd data
  hlon=zeros(192,94) ; hlat=zeros(192,94);
  for jj=1:94
    for ii=1:192
      hlon(ii,jj)=0+(ii-1)*1.875;
      hlat(ii,jj)=-88.542+(jj-1)*1.914893617;
      if(hlon(ii,jj) < 180 & hlon(ii,jj) > 0)
         uwnd(ii,jj)=NaN; vwnd(ii,jj)=NaN;
      end
    end
  end
  ind=find(hlon >= 180) ; hlon(ind)=hlon(ind)-360; clear ind;

  %% maping ice conc etc...
  m_proj('lambert','lat',[65 81],'lon',[-180 -120]); %ランベルト図法に指定
  hold on;

  [cs,hs]=m_contourf(lon,lat,ice,[10:1:100]);
  a=flipud(cptcmap('GMT_white-blue', 'mapping', 'showall','ncol',18));
  colormap(a) ; h=colorbar; caxis([10 100]);
  set(get(h,'ylabel'),'String','Sea ice concentration(%)','fontsize',12);
  set(hs,'LineStyle','none')

  %% landmask and coast line
  m_usercoast('cst','patch',[.4 .4 .4]); % 海岸線の描画

  %% Contour --- IBCAO
  [c2,h2]=m_contour(iblon,iblat,ibd2,[600:600:1800],'k','linewidth',1,'linest','-'); %海底地形のコンターライン

  %% Contour --- SLP
  [c,h]=m_contour(plon,plat,slp,[1008:2:1030],'k-','linewidth',2); %高気圧側
  clabel(c,h,'fontsize',12); h1 = clabel(c,h);
  set(h1,'FontName','Arial');
  set(h1,'Fontsize',12);
  [c,h]=m_contour(plon,plat,slp,[980:2:1006],'k-','linewidth',2,'linest','-.'); %低気圧側
  clabel(c,h,'fontsize',12); h1 = clabel(c,h);
  set(h1,'FontName','Arial');
  set(h1,'Fontsize',12);

  %% Vector mapping
  s=60 ; % ベクトルのスケールファクター
  ik=2 ; jk=2; %ベクトルデータの間引き加減
  lon2=slon(1:ik:159,1:jk:151);
  lat2=slat(1:ik:159,1:jk:151);
  ve2=mve(1:ik:159,1:jk:151);
  vn2=mvn(1:ik:159,1:jk:151);

  hh=m_quiver(lon2,lat2,ve2,vn2,4,'linewidth',1.)); % 氷速ベクトル描画
  set(hh,'Color',[0 0.6 1]);

  s2=30;
  ik=2 ; jk=2;
  hlon3=imresize(hlon,[192,188],'bilinear'); hlat3=imresize(hlat,[192,188],'bilinear');
  uwnd3=imresize(uwnd,[192,188],'bilinear'); vwnd3=imresize(vwnd,[192,188],'bilinear');

  hlon4=hlon3(1:ik:192,1:jk:188);
  hlat4=hlat3(1:ik:192,1:jk:188);
  uwnd4=uwnd3(1:ik:192,1:jk:188);
  vwnd4=vwnd3(1:ik:192,1:jk:188);

  m_vec(s2, -150, 67.7, 10, 0, 'r', 'key', 'Wind 10 m s^{-1}','EdgeColor','k'); % 風速ベクトルのリファレンス
  m_vec(s2,hlon4,hlat4,uwnd4,vwnd4,'r','shaftwidth',1.2,'headangle',45); % 風速ベクトル描画

  m_grid('box','fancy','tickdir','in');
  xlabel('Ice(shaded)+Uice(cyan)+SLP(blk)+WND(red)','fontsize',14);
  title([yr2 ' Frz-up season'],'fontsize',30,'fontweight','bold');

  hold off;
  print('-djpeg','-r200',ofl);
  close all;
end