script demo: map

Friday, May 3, 2013

[Matlab] remove the water bodies in a m_map figure

function m_nolakes(lakeColor)
% delete the water bodies patches in figures created by m_map
%           default lakeColor is white ([1 1 1])
% usage:
%           m_nolakes(lakeColor)
%           http://scriptdemo.blogspot.com

if nargin==0
   lakeColor=[1 1 1];
elseif nargin~=1
   help m_nolakes
   return
end

gshhstypes={'c','l','i','h','f'};
for nt=1:numel(gshhstypes)
    eval(['hw=findobj(gcf,''type'',''patch'',''tag'',''m_gshhs_',gshhstypes{nt},''',''facecolor'',lakeColor);']);
    if ~isempty(hw)
       delete(hw);
    end
end
hw=findobj(gcf,'type','patch','tag','m_coast','facecolor',lakeColor);
if ~isempty(hw)
   delete(hw);
end

Friday, July 23, 2010

matlab 3d plot with a map background

load korea
load geoid
t=0:pi/100:10*pi;
xx=100+50*sin(t);
yy=60+50*cos(t);
pcolor(map);shading interp; colorbar
colormap(demcmap(map))
cmap=colormap;
for N=1:length(xx)
C(N,1:3)=cmap(max(1,mod(N,64)),:);
end

hold on;
scatter3(xx(:),yy(:),t(:),10,C);
view(3)

Tuesday, March 23, 2010

matlab: shaded world map based on coast data

function FillCoast(Is360,FillColor)
% To show a shaded worldmap based on the coast data coming with matlab
% Not perfectly, but does work some way.
% Usage:
%       FillCoast(Is360,FillColor)
% Copyright: http://scriptdemo.blogspot.com
if (nargin<1)
     Is360=0;
     FillColor='r';
elseif   (nargin==1)
     if (isnumeric(Is360) && length(Is360)==1)
        FillColor='r';
     else
        FillColor=Is360;
        Is360=0;
     end
elseif   (nargin>2)
     help FillCoast
     return
end
load coast;
if   ~ishold  hold on; end
if   Is360
     long(long<0)=long(long<0)+360;
     lat1=lat;
     lat1(long<0.5)=nan;
     plot(long,lat1,'k'); clear lat1
end
IndNaN=find(isnan(long));
Ns=1;
for   N=1:length(IndNaN)-1
      Ne=IndNaN(N)-1;
      lonseg=long(Ns:Ne);
      latseg=lat(Ns:Ne);
      Ns=Ne+2;
     if (length(lonseg)>3)
        N360=find(lonseg>340);
        NZero=find(lonseg<10);
        switch ~isempty(N360) + ~isempty(NZero)
           case {0,1}
                if (~Is360 && isempty(find(latseg>-60)))
                   % Antarctic
                   lonseg=[lonseg; lonseg(end); lonseg(1)];
                   latseg=[latseg;         -90;       -90];
                end
                hf=patch(lonseg,latseg,FillColor); set(hf,'linestyle','none')
           case {2}
                if isempty(find(latseg>-60))
                % Antarctic
                  Nmin=find(lonseg==min(lonseg));
                  lonseg=[lonseg(1:Nmin-1); lonseg(Nmin-1); lonseg(Nmin); lonseg(Nmin:end)];
                  latseg=[latseg(1:Nmin-1);         -90;       -90; latseg(Nmin:end)];
                  hf=patch(lonseg,latseg,FillColor); set(hf,'linestyle','none')
               else
                 lonseg1=lonseg(lonseg>200 & lonseg<=360);
                 if ~isempty(lonseg1)
                    lonseg1=[lonseg1;lonseg1(1)];
                    latseg1=latseg(lonseg>200 & lonseg<=360);
                    latseg1=[latseg1;latseg1(1)];
                    hf=patch(lonseg1,latseg1,'b'); set(hf,'linestyle','none')
                 end;
                 lonseg2=lonseg(lonseg<195);
                 if ~isempty(lonseg2)
                    lonseg2=[lonseg2;lonseg2(1)];
                    latseg2=latseg(lonseg<195);  latseg2=[latseg2;latseg2(1)];
                    hf=patch(lonseg2,latseg2,'g'); set(hf,'linestyle','none')
                 end
               end
           otherwise
             disp('You should not see this info...STOP')
        end
     end
end