clear; close; clc

% latlim = [25 50]; lonlim = [-105 -67];  % Breiten- und Laengengrad-Bereich definieren
                                        % Teilbereich klappt
% latlim = [-90 90]; lonlim = [-180 180];  gesamte Darstellung haut nicht
% hin, 

% axe = usamap(latlim, lonlim);
% 
% landareas = shaperead('landareas.shp','UseGeoCoords', true, 'BoundingBox', [lonlim', latlim']);
% geoshow(landareas);

axe = worldmap('Europe');         % worldmap('World') geht für weltweit auch
% setm(axe, 'Origin', [0 180 0])  % benutzt mit worldmap('World')

landareas = shaperead('landareas', 'UseGeoCoords', true);
geoshow(axe, landareas, 'FaceColor', [0.5 0.7 0.5])

cities = shaperead('worldcities', 'UseGeoCoords', true);  % europecities gibt's nicht

x1 = (13.405);   % Berlin Longitude
y1 = (52.520);   % Berlin Latitude

x2 = (9.994);    % Hamburg Longitude
y2 = (53.551);   % Hamburg Latitude

%  Koordinaten von: https://www.latlong.net/place/hamburg-germany-8766.html
%  passt nicht ganz mit dem, was bei worldcities hinterlegt ist, also Pkt.
%  werden nicht genau getroffen

plotm(y1,x1, 'x', 'MarkerSize', 10, 'LineWidth', 2);  % hier muss zuerst y dann x
hold on;
plotm(y2,x2, 'x', 'MarkerSize', 10, 'LineWidth', 2);
textm(y1,x1, 'Berlin')
textm(y2,x2, 'Hamburg')

geoshow(cities, 'Marker', '.', 'Color', 'red')