OK, a full GMT example
# Extract a sub-region from a global grid. Using a low resolution to keep example light
Ggeog = grdcut("@earth_relief_10m", region=(-18,15,32,56));
convert to Lambert
Glamb = grdproject(Ggeog,proj="+proj=lcc+lat_1=39+lat_2=50+lat_0=45+lon_0=0");
Show them. Rafael, you are right about the extrapolation but the trick is to not do it and instead fill those nodes with NaNs.
imshow(Ggeog, shade=true, cmap=:geo, name="geogs.png")
imshow(Glamb, shade=true, cmap=:geo, name="lcc.png")
