    function WMTS_example

    %Open a new figure on monitor 2:
    HWIN=figure('pos',[2073 553 560 420]);

    %List WMTS layers that are available out of the box:
    WMTS_services=modelit.mapviewer.TileSourceFactory.allTileSources;

    disp('Available WMTS services')
    disp('======================================')
    disp(char(WMTS_services.displayName))
    disp('======================================')

    %In the example below we select 1 of the above. Feel free to experiment
    %with other WMTS services
    tilesource4example = modelit.mapviewer.TileSourceFactory.create(struct('displayName','Luchtfotos actueel (HR)'));
    % tilesource4example = modelit.mapviewer.TileSourceFactory.create(struct('displayName','USGSTopo'));
    % tilesource4example = modelit.mapviewer.TileSourceFactory.create(struct('displayName','OpenStreetMap'));

    m = modelit.mapviewer.MapViewer(HWIN,'tilesource',tilesource4example);

    %Center the map in the figure
    set(m.handle,'units','normalized','pos',[.1 .1 .8 .8]);

    %A Matlab figure with a zoomable map is now visible.
    %Its content is served by (for example) http://a.tile.openstreetmap.org
    %Continue by plotting some data of our own:

    %load linedata for the Dutch Motorway&Trunk roads network. These data are
    %reported in RD coordinates:
    data=load('highways.mat','-mat');

    if isa(tilesource4example,'modelit.mapviewer.WGS84TileSource')
        %tilesource4example is a WGS service. Conver RD to this coordinates
        %before plotting.
        xy=RD2WGS([data.x,data.y]);
        data.x=xy(:,1);
        data.y=xy(:,2);
    end

    %Plot the data
    h=line(data.x,data.y,'parent',m,'color','r','linew',2);

    %Select the Area to display (in this case: range of displayed line)
    setVisibleRectangle(m, [ min(data.x),min(data.y),max(data.x)-min(data.x),max(data.y)-min(data.y)]);

    %h is a native Matlab line object
    get(h)
    set(h,'color','b')
    set(h,'buttondown',@(o,e)disp('Button Pressed'))