how to plot country borders in regional data

I have a netCDF file in which I have a variable to plot at 12 Km grid spacing. The domain is regional and contains 3 countries i.e. Norway, Sweden, and Finland.

Is there a way to draw country borders?

Example plot is attached.

LST_JJActl  = ncread('/Volumes/YKMac/EnsAvg_NSF12KP3CTL_JJA_DTM.nc','LST');
pft=ncread('/Volumes/YKMac/NetCDF_Inputdata_P3/NSF12KP3CTL_CLM45_surface.nc','pft_2d');
nt=pft(:,:,3); bt=pft(:,:,9);
NETx=nan(size(LST_JJActl,1),size(LST_JJActl,2)); BDTx=nan(size(LST_JJActl,1),size(LST_JJActl,2));
NETx(1:end,1:end)=nt(2:end-2,2:end-2); BDTx(1:end,1:end)=bt(2:end-2,2:end-2);
NET=NETx*0.01; BDT=BDTx*0.01;
PFTtot=(NET+BDT); % Sum of NET and BDT pft fraction to apply on variables of grid level 

x0=10; y0=10; width=1200; height=600;

figure(1); 
set(gcf,'units','points','position',[x0,y0,width,height]);
%$*********************** PFT fraction ***********************$%

subplot(1,2,1);
h11=imagesc(rot90(NETx));axis on;grid; set(h11,'alphadata',~isnan(rot90(NETx)));
set(gca,'color', [0.5 0.5 0.5]); caxis([0 100]); colormap(bluewhitered(256));
yticklabels({'69','67','65','63','61','59','57','55'});
xticklabels({'0','5','10','15','20','25','30','35'});
ylabel('Lat');xlabel('Lon'); title('NET','fontsize',14);

subplot(1,2,2);
h12=imagesc(rot90(BDTx));axis on;grid; set(h12,'alphadata',~isnan(rot90(BDTx)));
set(gca,'color', [0.5 0.5 0.5]); caxis([0 100]); colormap(bluewhitered(256));
yticklabels({'69','67','65','63','61','59','57','55'});
xticklabels({'0','5','10','15','20','25','30','35'});
ylabel('Lat');xlabel('Lon'); title('BDT','fontsize',14);

hp1 = get(subplot(1,2,2),'Position');
hc1 = colorbar('Location','Manual','Position', [hp1(1)-0.062  hp1(2)+0.024  hp1(3)-0.32  hp1(4)-0.06]);
title(hc1,'%');
sgtitle('Fraction of PFT in RegCM+CLM, NET and BDT','fontweight','bold','fontsize',16); 
set(gcf,'inverthardcopy',false);
saveas(gcf,'/Volumes/YKMac/Figures_P3/FigP3_SI4_PFTfract.png');

Thanks.

Answer 1

[KSSV]

You may download your required countries boundary shape file from this link: https://gadm.org/download_country.html

You can read the data using shaperead and then plot.

[Yogesh Kumkar]

Thanks, I already have gdam40_FIN_shp, gdam40_NOR_shp, and gdam40_SWE_shp.

These folders are at the location that is in path. I get the error.

addpath('/Volumes/YKMac/MatlabFunctions/');
lat  = ncread('/Volumes/YKMac/EnsAvg_NSF12KP3CTL_JJA_DTM.nc','lat');
lon  = ncread('/Volumes/YKMac/EnsAvg_NSF12KP3CTL_JJA_DTM.nc','lon');
Norway=shaperead('Norway.shp');

Error using openShapeFiles>checkSHP (line 67)

Unable to open both Norway.shp and Norway.SHP.

Error in openShapeFiles (line 18)

[basename, ext] = checkSHP(basename,shapeExtensionProvided);

Error in shaperead (line 210)

= openShapeFiles(filename,'shaperead');

I am very new to this functionality.

All I want to do is plot a variable from the attached file and show country borders.

Thanks.

[KSSV]

Shape files will have few other supporting files, are these files too present in the same folder?

[Yogesh Kumkar]

yes. I have many files in each dir. e.g. cpg, dbf, prj, shp, shx,

Is there anything needed?

[KSSV]

Attach yoour shape files and complete code which you have tried.

[Yogesh Kumkar]

I tried to at least border Norway, but still get error.

Temp  = ncread('/Volumes/YogeshKumkarMac/PAPER_3/Data_P3/NetCDF_Inputdata_P3/TemperatureNSF.nc','ts');
lat  = ncread('/Volumes/YogeshKumkarMac/PAPER_3/Data_P3/NetCDF_Inputdata_P3/TemperatureNSF.nc','xlat');
lon  = ncread('/Volumes/YogeshKumkarMac/PAPER_3/Data_P3/NetCDF_Inputdata_P3/TemperatureNSF.nc','xlon');
figure(3);
h11=imagesc(rot90(Temp));axis on;grid; set(h11,'alphadata',~isnan(rot90(Temp)));
set(gca,'color', [0.5 0.5 0.5]); colormap(bluewhitered(256)); colorbar;
title('Temperature','fontsize',14);

Norway=shaperead('Norway.shp'); % does not work. Have addpath the location.
geoshow(lat,lon,'Color','black') % Also have iy and jx
mymap=pcolor(lon,lat,NETx');
mymap.EdgeAlpha=0;

set(gcf,'inverthardcopy',false);
saveas(gcf,'/Volumes/YogeshKumkarMac/PAPER_3/Figures_P3/Temperature.png');

I have skipped gadm40_NOR_2.shp and gadm40_NOR_0.shp because of size upload limit.

[KSSV]

It looks like you are giving wrong name for the shape file. Please check.

S = shaperead('gadm40_NOR_1.shp') ;
S = [[S.X]' [S.Y]'];
plot(S(:,1),S(:,2))

[Yogesh Kumkar]

Thanks, I could plot the borders of gdam40_NOR_0.shp, SWE, and FIN.

Looks good but, it does not overlay on my variable. May be something about projection and resolution I should do? The data I have is at 12 Km and at rotated mercator projection. Is there any way to exactly overlay the borders on variable?

N = shaperead('gadm40_NOR_0.shp') ;
S = shaperead('gadm40_SWE_0.shp') ;
F = shaperead('gadm40_FIN_0.shp') ;
figure(4);
plot(N.X, N.Y','k'); hold on;
plot(S.X, S.Y','k'); hold on;
plot(F.X, F.Y','k'); hold on;
title('Domain: Norway, Sweden, Finland','fontsize',14); 
set(gcf,'inverthardcopy',false);
saveas(gcf,'/Volumes/YogeshKumkarMac/PAPER_3/Figures_P3/DomainNSF.png');